{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Modeling and Simulation in Python\n",
    "\n",
    "Experiments with different ODE solvers\n",
    "\n",
    "Copyright 2019 Allen Downey\n",
    "\n",
    "License: [Creative Commons Attribution 4.0 International](https://creativecommons.org/licenses/by/4.0)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Configure Jupyter so figures appear in the notebook\n",
    "%matplotlib inline\n",
    "\n",
    "# Configure Jupyter to display the assigned value after an assignment\n",
    "%config InteractiveShell.ast_node_interactivity='last_expr_or_assign'\n",
    "\n",
    "# import functions from the modsim.py module\n",
    "from modsim import *"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>values</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>init</th>\n",
       "      <td>y    2\n",
       "dtype: int64</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>t_0</th>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>t_end</th>\n",
       "      <td>3</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "init     y    2\n",
       "dtype: int64\n",
       "t_0                        1\n",
       "t_end                      3\n",
       "dtype: object"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "init = State(y = 2)\n",
    "system = System(init=init, t_0=1, t_end=3)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [],
   "source": [
    "def slope_func(state, t, system):\n",
    "    [y] = state\n",
    "    dydt = y + t\n",
    "    return [dydt]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [],
   "source": [
    "results, details = run_euler(system, slope_func)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "24.973714732487576"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "get_last_value(results.y)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Glucose minimal model\n",
    "\n",
    "Read the data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "data = pd.read_csv('data/glucose_insulin.csv', index_col='time');"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Interpolate the insulin data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<scipy.interpolate.interpolate.interp1d at 0x7f0b8f890908>"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "I = interpolate(data.insulin)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Initialize the parameters"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "1e-05"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "G0 = 290\n",
    "k1 = 0.03\n",
    "k2 = 0.02\n",
    "k3 = 1e-05"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To estimate basal levels, we'll use the concentrations at `t=0`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "11"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Gb = data.glucose[0]\n",
    "Ib = data.insulin[0]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Create the initial condtions."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>values</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>G</th>\n",
       "      <td>290</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>X</th>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "G    290\n",
       "X      0\n",
       "dtype: int64"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "init = State(G=G0, X=0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Make the `System` object."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "182"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "t_0 = get_first_label(data)\n",
    "t_end = get_last_label(data)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>values</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>G0</th>\n",
       "      <td>290</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>k1</th>\n",
       "      <td>0.03</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>k2</th>\n",
       "      <td>0.02</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>k3</th>\n",
       "      <td>1e-05</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>init</th>\n",
       "      <td>G    290\n",
       "X      0\n",
       "dtype: int64</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>Gb</th>\n",
       "      <td>92</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>Ib</th>\n",
       "      <td>11</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>I</th>\n",
       "      <td>&lt;scipy.interpolate.interpolate.interp1d object...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>t_0</th>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>t_end</th>\n",
       "      <td>182</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>dt</th>\n",
       "      <td>2</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "G0                                                     290\n",
       "k1                                                    0.03\n",
       "k2                                                    0.02\n",
       "k3                                                   1e-05\n",
       "init                        G    290\n",
       "X      0\n",
       "dtype: int64\n",
       "Gb                                                      92\n",
       "Ib                                                      11\n",
       "I        <scipy.interpolate.interpolate.interp1d object...\n",
       "t_0                                                      0\n",
       "t_end                                                  182\n",
       "dt                                                       2\n",
       "dtype: object"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "system = System(G0=G0, k1=k1, k2=k2, k3=k3,\n",
    "                init=init, Gb=Gb, Ib=Ib, I=I,\n",
    "                t_0=t_0, t_end=t_end, dt=2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "def update_func(state, t, system):\n",
    "    \"\"\"Updates the glucose minimal model.\n",
    "    \n",
    "    state: State object\n",
    "    t: time in min\n",
    "    system: System object\n",
    "    \n",
    "    returns: State object\n",
    "    \"\"\"\n",
    "    G, X = state\n",
    "    k1, k2, k3 = system.k1, system.k2, system.k3 \n",
    "    I, Ib, Gb = system.I, system.Ib, system.Gb\n",
    "    dt = system.dt\n",
    "    \n",
    "    dGdt = -k1 * (G - Gb) - X*G\n",
    "    dXdt = k3 * (I(t) - Ib) - k2 * X\n",
    "    \n",
    "    G += dGdt * dt\n",
    "    X += dXdt * dt\n",
    "\n",
    "    return State(G=G, X=X)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "def run_simulation(system, update_func):\n",
    "    \"\"\"Runs a simulation of the system.\n",
    "        \n",
    "    system: System object\n",
    "    update_func: function that updates state\n",
    "    \n",
    "    returns: TimeFrame\n",
    "    \"\"\"\n",
    "    t_0, t_end, dt = system.t_0, system.t_end, system.dt\n",
    "    \n",
    "    frame = TimeFrame(columns=init.index)\n",
    "    frame.row[t_0] = init\n",
    "    ts = linrange(t_0, t_end, dt)\n",
    "    \n",
    "    for t in ts:\n",
    "        frame.row[t+dt] = update_func(frame.row[t], t, system)\n",
    "    \n",
    "    return frame"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "CPU times: user 166 ms, sys: 4.08 ms, total: 170 ms\n",
      "Wall time: 167 ms\n"
     ]
    }
   ],
   "source": [
    "%time results = run_simulation(system, update_func);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>G</th>\n",
       "      <th>X</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>290</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>278.12</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>266.953</td>\n",
       "      <td>0.0003</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6</th>\n",
       "      <td>256.295</td>\n",
       "      <td>0.002668</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8</th>\n",
       "      <td>245.07</td>\n",
       "      <td>0.00404128</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>10</th>\n",
       "      <td>233.905</td>\n",
       "      <td>0.00467963</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>12</th>\n",
       "      <td>223.202</td>\n",
       "      <td>0.00525244</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>14</th>\n",
       "      <td>212.985</td>\n",
       "      <td>0.00572235</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>16</th>\n",
       "      <td>203.288</td>\n",
       "      <td>0.00609345</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>18</th>\n",
       "      <td>194.133</td>\n",
       "      <td>0.00632971</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>20</th>\n",
       "      <td>185.548</td>\n",
       "      <td>0.00648986</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>22</th>\n",
       "      <td>177.527</td>\n",
       "      <td>0.00661026</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>24</th>\n",
       "      <td>170.048</td>\n",
       "      <td>0.00672585</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>26</th>\n",
       "      <td>163.078</td>\n",
       "      <td>0.00681282</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>28</th>\n",
       "      <td>156.591</td>\n",
       "      <td>0.00687231</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>30</th>\n",
       "      <td>150.563</td>\n",
       "      <td>0.00692941</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>32</th>\n",
       "      <td>144.963</td>\n",
       "      <td>0.00700824</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>34</th>\n",
       "      <td>139.753</td>\n",
       "      <td>0.00710791</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>36</th>\n",
       "      <td>134.901</td>\n",
       "      <td>0.00717159</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>38</th>\n",
       "      <td>130.392</td>\n",
       "      <td>0.00720073</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>40</th>\n",
       "      <td>126.211</td>\n",
       "      <td>0.0071967</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>42</th>\n",
       "      <td>122.342</td>\n",
       "      <td>0.00716083</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>44</th>\n",
       "      <td>118.769</td>\n",
       "      <td>0.0070944</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>46</th>\n",
       "      <td>115.478</td>\n",
       "      <td>0.00700262</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>48</th>\n",
       "      <td>112.452</td>\n",
       "      <td>0.00688652</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>50</th>\n",
       "      <td>109.676</td>\n",
       "      <td>0.00674706</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>52</th>\n",
       "      <td>107.135</td>\n",
       "      <td>0.00658517</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>54</th>\n",
       "      <td>104.816</td>\n",
       "      <td>0.00640177</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>56</th>\n",
       "      <td>102.705</td>\n",
       "      <td>0.0062257</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>58</th>\n",
       "      <td>100.784</td>\n",
       "      <td>0.00605667</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>...</th>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>124</th>\n",
       "      <td>86.3907</td>\n",
       "      <td>0.00109474</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>126</th>\n",
       "      <td>86.5381</td>\n",
       "      <td>0.000972951</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>128</th>\n",
       "      <td>86.6974</td>\n",
       "      <td>0.000858033</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>130</th>\n",
       "      <td>86.8668</td>\n",
       "      <td>0.000749712</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>132</th>\n",
       "      <td>87.0445</td>\n",
       "      <td>0.000647723</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>134</th>\n",
       "      <td>87.2291</td>\n",
       "      <td>0.000551815</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>136</th>\n",
       "      <td>87.4191</td>\n",
       "      <td>0.000461742</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>138</th>\n",
       "      <td>87.6132</td>\n",
       "      <td>0.000377272</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>140</th>\n",
       "      <td>87.8103</td>\n",
       "      <td>0.000298181</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>142</th>\n",
       "      <td>88.0093</td>\n",
       "      <td>0.000224254</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>144</th>\n",
       "      <td>88.2093</td>\n",
       "      <td>0.000155284</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>146</th>\n",
       "      <td>88.4093</td>\n",
       "      <td>8.90726e-05</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>148</th>\n",
       "      <td>88.609</td>\n",
       "      <td>2.55097e-05</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>150</th>\n",
       "      <td>88.808</td>\n",
       "      <td>-3.55107e-05</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>152</th>\n",
       "      <td>89.0058</td>\n",
       "      <td>-9.40903e-05</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>154</th>\n",
       "      <td>89.2022</td>\n",
       "      <td>-0.000150327</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>156</th>\n",
       "      <td>89.3969</td>\n",
       "      <td>-0.000204314</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>158</th>\n",
       "      <td>89.5896</td>\n",
       "      <td>-0.000256141</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>160</th>\n",
       "      <td>89.7801</td>\n",
       "      <td>-0.000305895</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>162</th>\n",
       "      <td>89.9682</td>\n",
       "      <td>-0.00035366</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>164</th>\n",
       "      <td>90.1538</td>\n",
       "      <td>-0.000399513</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>166</th>\n",
       "      <td>90.3366</td>\n",
       "      <td>-0.000445533</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>168</th>\n",
       "      <td>90.5169</td>\n",
       "      <td>-0.000491711</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>170</th>\n",
       "      <td>90.6949</td>\n",
       "      <td>-0.000538043</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>172</th>\n",
       "      <td>90.8708</td>\n",
       "      <td>-0.000584521</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>174</th>\n",
       "      <td>91.0448</td>\n",
       "      <td>-0.00063114</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>176</th>\n",
       "      <td>91.217</td>\n",
       "      <td>-0.000677895</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>178</th>\n",
       "      <td>91.3877</td>\n",
       "      <td>-0.000724779</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>180</th>\n",
       "      <td>91.5569</td>\n",
       "      <td>-0.000771788</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>182</th>\n",
       "      <td>91.7248</td>\n",
       "      <td>-0.000818916</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>92 rows × 2 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "           G            X\n",
       "0        290            0\n",
       "2     278.12            0\n",
       "4    266.953       0.0003\n",
       "6    256.295     0.002668\n",
       "8     245.07   0.00404128\n",
       "10   233.905   0.00467963\n",
       "12   223.202   0.00525244\n",
       "14   212.985   0.00572235\n",
       "16   203.288   0.00609345\n",
       "18   194.133   0.00632971\n",
       "20   185.548   0.00648986\n",
       "22   177.527   0.00661026\n",
       "24   170.048   0.00672585\n",
       "26   163.078   0.00681282\n",
       "28   156.591   0.00687231\n",
       "30   150.563   0.00692941\n",
       "32   144.963   0.00700824\n",
       "34   139.753   0.00710791\n",
       "36   134.901   0.00717159\n",
       "38   130.392   0.00720073\n",
       "40   126.211    0.0071967\n",
       "42   122.342   0.00716083\n",
       "44   118.769    0.0070944\n",
       "46   115.478   0.00700262\n",
       "48   112.452   0.00688652\n",
       "50   109.676   0.00674706\n",
       "52   107.135   0.00658517\n",
       "54   104.816   0.00640177\n",
       "56   102.705    0.0062257\n",
       "58   100.784   0.00605667\n",
       "..       ...          ...\n",
       "124  86.3907   0.00109474\n",
       "126  86.5381  0.000972951\n",
       "128  86.6974  0.000858033\n",
       "130  86.8668  0.000749712\n",
       "132  87.0445  0.000647723\n",
       "134  87.2291  0.000551815\n",
       "136  87.4191  0.000461742\n",
       "138  87.6132  0.000377272\n",
       "140  87.8103  0.000298181\n",
       "142  88.0093  0.000224254\n",
       "144  88.2093  0.000155284\n",
       "146  88.4093  8.90726e-05\n",
       "148   88.609  2.55097e-05\n",
       "150   88.808 -3.55107e-05\n",
       "152  89.0058 -9.40903e-05\n",
       "154  89.2022 -0.000150327\n",
       "156  89.3969 -0.000204314\n",
       "158  89.5896 -0.000256141\n",
       "160  89.7801 -0.000305895\n",
       "162  89.9682  -0.00035366\n",
       "164  90.1538 -0.000399513\n",
       "166  90.3366 -0.000445533\n",
       "168  90.5169 -0.000491711\n",
       "170  90.6949 -0.000538043\n",
       "172  90.8708 -0.000584521\n",
       "174  91.0448  -0.00063114\n",
       "176   91.217 -0.000677895\n",
       "178  91.3877 -0.000724779\n",
       "180  91.5569 -0.000771788\n",
       "182  91.7248 -0.000818916\n",
       "\n",
       "[92 rows x 2 columns]"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "results"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Numerical solution\n",
    "\n",
    "In the previous chapter, we approximated the differential equations with difference equations, and solved them using `run_simulation`.\n",
    "\n",
    "In this chapter, we solve the differential equation numerically using `run_euler`...\n",
    "\n",
    "Instead of an update function, we provide a slope function that evaluates the right-hand side of the differential equations.  We don't have to do the update part; the solver does it for us."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [],
   "source": [
    "def slope_func(state, t, system):\n",
    "    \"\"\"Computes derivatives of the glucose minimal model.\n",
    "    \n",
    "    state: State object\n",
    "    t: time in min\n",
    "    system: System object\n",
    "    \n",
    "    returns: derivatives of G and X\n",
    "    \"\"\"\n",
    "    G, X = state\n",
    "    k1, k2, k3 = system.k1, system.k2, system.k3 \n",
    "    I, Ib, Gb = system.I, system.Ib, system.Gb\n",
    "    \n",
    "    dGdt = -k1 * (G - Gb) - X*G\n",
    "    dXdt = k3 * (I(t) - Ib) - k2 * X\n",
    "    \n",
    "    return dGdt, dXdt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can test the slope function with the initial conditions."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(-5.9399999999999995, 0.0)"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "slope_func(init, 0, system)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here's how we run the ODE solver."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "CPU times: user 276 ms, sys: 368 µs, total: 277 ms\n",
      "Wall time: 274 ms\n"
     ]
    }
   ],
   "source": [
    "system = System(G0=G0, k1=k1, k2=k2, k3=k3,\n",
    "                init=init, Gb=Gb, Ib=Ib, I=I,\n",
    "                t_0=t_0, t_end=t_end, dt=1)\n",
    "\n",
    "%time results2, details = run_euler(system, slope_func)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`results` is a `TimeFrame` with one row for each time step and one column for each state variable:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>G</th>\n",
       "      <th>X</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>290</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>284.06</td>\n",
       "      <td>0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>278.298</td>\n",
       "      <td>7.5e-05</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>272.688</td>\n",
       "      <td>0.0002235</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>267.207</td>\n",
       "      <td>0.00088903</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5</th>\n",
       "      <td>261.713</td>\n",
       "      <td>0.00206125</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6</th>\n",
       "      <td>256.082</td>\n",
       "      <td>0.00298502</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7</th>\n",
       "      <td>250.395</td>\n",
       "      <td>0.00366532</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8</th>\n",
       "      <td>244.726</td>\n",
       "      <td>0.00416202</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9</th>\n",
       "      <td>239.125</td>\n",
       "      <td>0.00447878</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>10</th>\n",
       "      <td>233.641</td>\n",
       "      <td>0.0047792</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>11</th>\n",
       "      <td>228.275</td>\n",
       "      <td>0.00506362</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>12</th>\n",
       "      <td>223.031</td>\n",
       "      <td>0.00532235</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>13</th>\n",
       "      <td>217.913</td>\n",
       "      <td>0.0055559</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>14</th>\n",
       "      <td>212.925</td>\n",
       "      <td>0.00576478</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>15</th>\n",
       "      <td>208.069</td>\n",
       "      <td>0.00594948</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>16</th>\n",
       "      <td>203.349</td>\n",
       "      <td>0.0061005</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>17</th>\n",
       "      <td>198.768</td>\n",
       "      <td>0.00621849</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>18</th>\n",
       "      <td>194.329</td>\n",
       "      <td>0.00631745</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>19</th>\n",
       "      <td>190.032</td>\n",
       "      <td>0.00639777</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>20</th>\n",
       "      <td>185.875</td>\n",
       "      <td>0.00645981</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>21</th>\n",
       "      <td>181.858</td>\n",
       "      <td>0.00652061</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>22</th>\n",
       "      <td>177.976</td>\n",
       "      <td>0.0065802</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>23</th>\n",
       "      <td>174.226</td>\n",
       "      <td>0.0066386</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>24</th>\n",
       "      <td>170.603</td>\n",
       "      <td>0.00668983</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>25</th>\n",
       "      <td>167.103</td>\n",
       "      <td>0.00673403</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>26</th>\n",
       "      <td>163.725</td>\n",
       "      <td>0.00677135</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>27</th>\n",
       "      <td>160.465</td>\n",
       "      <td>0.00680192</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>28</th>\n",
       "      <td>157.319</td>\n",
       "      <td>0.00682588</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>29</th>\n",
       "      <td>154.286</td>\n",
       "      <td>0.00685537</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>...</th>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>153</th>\n",
       "      <td>89.2245</td>\n",
       "      <td>-0.000116549</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>154</th>\n",
       "      <td>89.3181</td>\n",
       "      <td>-0.000144218</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>155</th>\n",
       "      <td>89.4115</td>\n",
       "      <td>-0.000171334</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>156</th>\n",
       "      <td>89.5045</td>\n",
       "      <td>-0.000197907</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>157</th>\n",
       "      <td>89.597</td>\n",
       "      <td>-0.000223949</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>158</th>\n",
       "      <td>89.6892</td>\n",
       "      <td>-0.00024947</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>159</th>\n",
       "      <td>89.7809</td>\n",
       "      <td>-0.000274481</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>160</th>\n",
       "      <td>89.8721</td>\n",
       "      <td>-0.000298991</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>161</th>\n",
       "      <td>89.9628</td>\n",
       "      <td>-0.000323011</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>162</th>\n",
       "      <td>90.053</td>\n",
       "      <td>-0.000346551</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>163</th>\n",
       "      <td>90.1426</td>\n",
       "      <td>-0.00036962</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>164</th>\n",
       "      <td>90.2316</td>\n",
       "      <td>-0.000392728</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>165</th>\n",
       "      <td>90.3201</td>\n",
       "      <td>-0.000415873</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>166</th>\n",
       "      <td>90.4081</td>\n",
       "      <td>-0.000439056</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>167</th>\n",
       "      <td>90.4955</td>\n",
       "      <td>-0.000462274</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>168</th>\n",
       "      <td>90.5825</td>\n",
       "      <td>-0.000485529</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>169</th>\n",
       "      <td>90.669</td>\n",
       "      <td>-0.000508818</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>170</th>\n",
       "      <td>90.7551</td>\n",
       "      <td>-0.000532142</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>171</th>\n",
       "      <td>90.8407</td>\n",
       "      <td>-0.000555499</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>172</th>\n",
       "      <td>90.926</td>\n",
       "      <td>-0.000578889</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>173</th>\n",
       "      <td>91.0108</td>\n",
       "      <td>-0.000602311</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>174</th>\n",
       "      <td>91.0953</td>\n",
       "      <td>-0.000625765</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>175</th>\n",
       "      <td>91.1795</td>\n",
       "      <td>-0.00064925</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>176</th>\n",
       "      <td>91.2633</td>\n",
       "      <td>-0.000672765</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>177</th>\n",
       "      <td>91.3468</td>\n",
       "      <td>-0.00069631</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>178</th>\n",
       "      <td>91.43</td>\n",
       "      <td>-0.000719883</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>179</th>\n",
       "      <td>91.5129</td>\n",
       "      <td>-0.000743486</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>180</th>\n",
       "      <td>91.5955</td>\n",
       "      <td>-0.000767116</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>181</th>\n",
       "      <td>91.6779</td>\n",
       "      <td>-0.000790774</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>182</th>\n",
       "      <td>91.7601</td>\n",
       "      <td>-0.000814458</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>183 rows × 2 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "           G            X\n",
       "0        290            0\n",
       "1     284.06            0\n",
       "2    278.298      7.5e-05\n",
       "3    272.688    0.0002235\n",
       "4    267.207   0.00088903\n",
       "5    261.713   0.00206125\n",
       "6    256.082   0.00298502\n",
       "7    250.395   0.00366532\n",
       "8    244.726   0.00416202\n",
       "9    239.125   0.00447878\n",
       "10   233.641    0.0047792\n",
       "11   228.275   0.00506362\n",
       "12   223.031   0.00532235\n",
       "13   217.913    0.0055559\n",
       "14   212.925   0.00576478\n",
       "15   208.069   0.00594948\n",
       "16   203.349    0.0061005\n",
       "17   198.768   0.00621849\n",
       "18   194.329   0.00631745\n",
       "19   190.032   0.00639777\n",
       "20   185.875   0.00645981\n",
       "21   181.858   0.00652061\n",
       "22   177.976    0.0065802\n",
       "23   174.226    0.0066386\n",
       "24   170.603   0.00668983\n",
       "25   167.103   0.00673403\n",
       "26   163.725   0.00677135\n",
       "27   160.465   0.00680192\n",
       "28   157.319   0.00682588\n",
       "29   154.286   0.00685537\n",
       "..       ...          ...\n",
       "153  89.2245 -0.000116549\n",
       "154  89.3181 -0.000144218\n",
       "155  89.4115 -0.000171334\n",
       "156  89.5045 -0.000197907\n",
       "157   89.597 -0.000223949\n",
       "158  89.6892  -0.00024947\n",
       "159  89.7809 -0.000274481\n",
       "160  89.8721 -0.000298991\n",
       "161  89.9628 -0.000323011\n",
       "162   90.053 -0.000346551\n",
       "163  90.1426  -0.00036962\n",
       "164  90.2316 -0.000392728\n",
       "165  90.3201 -0.000415873\n",
       "166  90.4081 -0.000439056\n",
       "167  90.4955 -0.000462274\n",
       "168  90.5825 -0.000485529\n",
       "169   90.669 -0.000508818\n",
       "170  90.7551 -0.000532142\n",
       "171  90.8407 -0.000555499\n",
       "172   90.926 -0.000578889\n",
       "173  91.0108 -0.000602311\n",
       "174  91.0953 -0.000625765\n",
       "175  91.1795  -0.00064925\n",
       "176  91.2633 -0.000672765\n",
       "177  91.3468  -0.00069631\n",
       "178    91.43 -0.000719883\n",
       "179  91.5129 -0.000743486\n",
       "180  91.5955 -0.000767116\n",
       "181  91.6779 -0.000790774\n",
       "182  91.7601 -0.000814458\n",
       "\n",
       "[183 rows x 2 columns]"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "results2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Plotting the results from `run_simulation` and `run_euler`, we can see that they are not very different."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[<matplotlib.lines.Line2D at 0x7f0b8f7954e0>]"
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXwAAAD9CAYAAAC/fMwDAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjAsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+17YcXAAAgAElEQVR4nO3deZCc9X3n8ffz9DF9zKUeaUbX6EBIPyOMwCDwgQ0E2RCyXid2bG9tIGsKvJt1eXdrw6Zq15Uqe51NTLzrpMImqeDsuioHELKJFjtgG2HkRVgG2wgkJCT4jSTQaCTNpek5u6en++l+9o/uGbUGHSPm6Jnuz6uqq/u5er79zMynf/17fv08ju/7iIhI9XMrXYCIiCwMBb6ISI1Q4IuI1AgFvohIjVDgi4jUiGClC7gYY0wdcDPQDeQrXI6IyFIRAFYBr1hrJ8oXLNrApxj2P6l0ESIiS9THgL3lMxZz4HcDPP7446xcubLStYiILAk9PT3ce++9UMrQcos58PMAK1euZO3atZWuRURkqXlXV7gO2oqI1AgFvohIjVDgi4jUiBn14RtjPgl8A9gI9AH/3Vr77dLQyVEgW7b6S9bau0rbfb603SpgD3C/tbZvDusXEZEZumzgG2NWAf8IfNpa+0NjzI3AT40xr1D8hJC01r5rGI0xZivwHeAeYB/wTeBJ4M45rP+SOjqT7N7XRW8yTVsixo7t7WxZn1ioHy8isqhctkvHWtsNrCiFvQu0AB7Flv1NwIGLbHof8LS1dq+1NgN8BbjVGLN5bkq/tI7OJI/tsoyksrQ0RRhJZXlsl6WjM7kQP15EZNGZUR++tXbUGBMDJoDngD+31h4FbgRajTEHjTG9xph/MMasKW22FThS9hxpoAu4bk5fwUXs3tdFPBIkHg3hOA7xaIh4JMjufV0L8eNFRBadKzlomwHiFL8B+4Ax5kEgBfwU2AEYYBx4qrR+PZCe9hxpIDabgmeqN5kmFjm/xyoWCdKbnF6SiEhtmPEXr6y1BYoHZ/cZY/4S+FVr7afK1zHGPAT0G2PaKb4ZRKc9TQwYm13JM9OWiDGSyhKPhqbmpTMebYkFeb8REVl0LtvCN8bcbox5ddrsOmDIGPN7xphryuaHS/cZit05pux5YsA6yrp55tOO7e2kMh6p8Ry+75Maz5HKeOzY3r4QP15EZNGZSQv/ALCm1Hp/BPgg8CDwaeC3ge3GmN8orfsI8H1rbb8x5glgrzHmDuBl4GFgv7W2Y45fwwVtWZ/gvrvNeaN0Pn37Jo3SEZGaddnAt9YOG2N+BfifwNcoHnj9orV2jzHmjdL8Y6Xn+j7wb0rbHTLGPAA8CqwBfg58bl5exUVsWZ9QwIuIlMyoD99a+xrw0QvMHwDuvcR2O4Gd77k6ERGZMzq1gohIjVDgi4jUCAW+iEiNUOCLiNQIBb6ISI1Q4IuI1AgFvohIjVDgi4jUCAW+iEiNUOCLiNQIBb6ISI1Q4IuI1AgFvohIjVDgi4jUCAW+iEiNUOCLiNQIBb6ISI1Q4IuI1AgFvohIjVDgi4jUCAW+iEiNUOCLiNQIBb6ISI1Q4IuI1AgFvohIjVDgi4jUiOBMVjLGfBL4BrAR6AP+u7X228aYMPBnwGeBPPDH1tqHy7b7fGm7VcAe4H5rbd/cvgQREZmJy7bwjTGrgH8E/rO1tgH4HPAnxpgbga8DBtgE3Ax8wRjzr0rbbQW+A9wPtABHgSfn4TWIiMgMXDbwrbXdwApr7Q+NMS7F8PaAUeALwB9YawettSeAbwG/Vdr0PuBpa+1ea20G+ApwqzFm8zy8DhERuYwZ9eFba0eNMTFgAngO+HOgn2JXzZGyVd8Cris93lq+zFqbBrrKlouIyAKaUR9+SQaIA9uAHwDjpfnpsnXSQKz0uH7asunLRURkAc048K21BSAL7DPG/CWwvbQoWrZaDBgrPU5NWzZ9uYiILKCZHLS93Rjz6rTZdcAg0EPxoO2k93GuG+dI+bJSl9A6zu8CEhGRBTKTFv4BYI0x5iHgEeCDwIPApykG/teMMQcpduH8TmkdgCeAvcaYO4CXgYeB/dbajjl9BSIiMiMzGaUzDPwK8BkgCfwl8EVr7R7gq8AbwGHgFWAn8Ghpu0PAA6Xps8C1FId0iohIBcyoD99a+xrw0QvMzwBfLt0utN1Oim8CC2qwvx+/4JNoa13oHy0ismhdySidJaH3ZBdvPvM4FApsuPPX2LB169Syjs4ku/d10ZtM05aIsWN7O1vWJypYrYjIwqm6c+mMDJyFQh7w6Xr1p1PzOzqTPLbLMpLK0tIUYSSV5bFdlo7OZOWKFRFZQFUX+Ks2bQKn+LLyw32cPdMNwO59XcQjQeLREI7jEI+GiEeC7N7XVclyRUQWTNUFfn1jI+G2jVPTJ17fB0BvMk0scn4PViwSpDc5/bthIiLVqeoCH2D1+2+YepzqfAsvm6UtESOd8c5bL53xaEvoi78iUhuqMvDbN2+BunoAfC/LicNvsGN7O6mMR2o8h+/7pMZzpDIeO7a3V7haEZGFUZWBHwi4NG16/9R071sH2LI+wX13GxrjYQaGMzTGw9x3t9EoHRGpGVU3LHPSxhtu4sCRnwM+uYFuhvv72LK+VQEvIjWrKlv4AM2JZgLLz3XXvHNg+umARERqS9UGPsDKrecO3o68c4S8511ibRGR6lbVgb9x6zX44eIonEI2Q+fhgxWuSESkcqo68IPBIA0bzx287Tm8v4LViIhUVlUHPsCmG7fjl755m032MHDmTIUrEhGpjKoP/GUtywi1bpia7jzwSuWKERGpoKoPfIC127ZPPR7tfItsZvwSa4uIVKeaCPz2qzfhR5sA8PMe7xxQX76I1J6aCPxAwCWx5dwQzf63DuD7fgUrEhFZeDUR+ACbb7wR3y1+sdgbG6TnnbcrXJGIyMKqmcCPxWNE1myZmu7cr4O3IlJbaibwATbeeMvU40zP24wODlWwGhGRhVVTgb+yfQ2B5pXFCd/n2Gtq5YtI7aipwAdYee2NU4+Hjh3Cy+UqWI2IyMKpucDfeN11OOEoAE5unOOHDle4IhGRhVFzgR8MBll29bnz63QfflVDNEWkJtRc4ANcfdMtOI5TnBju4czJ05UtSERkAdRk4Meamomt3jg1feLAvgpWIyKyMGZ0iUNjzCeAPwQ2A33A/7DWftsYUweMAtmy1V+y1t5V2u7zwDeAVcAe4H5rbd8c1v+e5Vdey/OvjzPshWkcGKLenOH6962udFkiIvPmsoFvjGkHdgJfAL4H3ATsMsacAAaApLV25QW22wp8B7gH2Ad8E3gSuHOOan/POjqT/ODgGBNuhIbABJm8y1997wBfikZ0zVsRqVoz6dLZADxhrX3KWluw1r4CvADcSjH8D1xku/uAp621e621GeArwK3GmM2zL3t2du/rIh4N0dzSguNAxC3gjg/xo1dOVro0EZF5c9kWvrX2J8BPJqeNMQngY8DfAr8MtBpjDgJtwIvAf7TWnga2UmzZTz5P2hjTBVwHHJ3LF3GlepNpWpoi+KEW0me7wfep8yfoPDVQybJERObVFR20NcY0Af8E/Jxi904K+CmwAzDAOPBUafV6ID3tKdJAbBb1zom2RIx0xsMNBgnXNwMw4bvUZQcrXJmIyPyZ0UFbAGPMFoohfwS411pbAB6ats5DQH+p3z8FRKc9TQwYm1XFc2DH9nYe22UBqG9ZzpnhYbKFADe5XQwlh2lONFW4QhGRuTejFr4x5jaKrfrvAp8t9cljjPk9Y8w1ZauGS/cZim8Mpuw5YsC60vyK2rI+wX13GxrjYUZzAeqjYT7S3MfqSJqjr2mIpohUp5mM0tkEPAP8rrX2T6ct3gZsN8b8Rmn6EeD71tp+Y8wTwF5jzB3Ay8DDwH5rbcecVT8LW9YnpkbkvHMgTOfe4vnxR98+RM67g1AwUMnyRETm3Ey6dL4MNAAPG2MeLpv/58CDwP8EjpWe6/vAvwGw1h4yxjwAPAqsofgJ4XNzV/rcyTS0s3twLUMTAZqCWfIvvMYnPn5zpcsSEZlTMxml8xDT+uqnufcS2+6kOIZ/0eroTPJ3u4/jRZbR4A0wXgjwT3s7Wb95k8bki0hVqclTK5Tbva+LeCRI66rWqTH54UKKp/csip4nEZE5U/OB35tME4sECdVFCEbrAahzCnR19Va4MhGRuVXzgT85Jh8gnlgBFMfkx71BxjMTlSxNRGRO1Xzg79jeTirjkRrPEW1sZMIJky0EeH88ybH9r1e6PBGROVPzgV8+Jn9gZIJlLQk+0tzHqsg4/faALo4iIlVjxt+0rWblY/InUqP87K8P4BfAHeun6+1O1m3aUNkCRUTmQM238KerizcQX3v11PTJ/T+vYDUiInNHgX8BG2784NTjXO/bDJ7VSdVEZOlT4F/A8rXrqEsUr+ni+AWOvapWvogsfQr8i1i17dypFcbeeYMJDdEUkSVOgX8R67ZeS6D0RSzHm+Do/tcqXJGIyOwo8C/CdQO0mA9MTZ998zUKhUIFKxIRmR0F/iVsumk7TqA0cjU9ROfiOLOziMh7osC/hLpojPj6c9d3OX3gFxWsRkRkdhT4l3H1zR8GHAC8gS56T3dXtiARkfdIgX8ZzStWEFm5fmr67VdermA1IiLvnU6tcAEdnUl27+uiN5mmLRHj+tXboOcEAJnTHYwOj9DQ1FjZIkVErpBa+NN0dCZ5bJdlJJWlpSnCSCrLs4dG6Q9MfhErj/3FzypcpYjIlVPgTzN5Bax4NITjOMSjIeLREMfdDVPrjB5/nYnxTOWKFBF5DxT400xeAWvS0GiGE90jHDo1we6hdXRnojhelo5XNWJHRJYWBf405VfAGhrNcOzUMJkJj8b6MIXYcl4aaqU7E2XgyKt4uVyFqxURmTkF/jTlV8A63Z/CAXxgbWsDbWvaCAfg8FgzZNMc23+g0uWKiMyYAn+a8itgjaSyROoCbG5vpqm+Dtd1aVi2jGEvDEDvGz8nn9fpFkRkadCwzAsovwLWSCpLPBqaWhZuStA03A+Anx7mncNvcPW2bRWpU0TkSqiFfwnl3Tu+79N9doyjp0YZcJbx/NlVdGeinNn/sk6qJiJLggL/Esq7dzq7RzjTn2LNijhbrlrJeCHIS0OtnO5PcfKtNytdqojIZc2oS8cY8wngD4HNQB/wP6y13zbGhIE/Az4L5IE/ttY+XLbd54FvAKuAPcD91tq+uX0J82uye+cvdr7OimXnunfqm5oYGx7i8Fgz7a/+lPXXbMVxnApXKyJycZdt4Rtj2oGdwO8DzcC/BB42xtwNfB0wwCbgZuALxph/VdpuK/Ad4H6gBTgKPDn3L2FhTB+f39y2krDjM+yF8Yb7OdVhK1idiMjlzaRLZwPwhLX2KWttwVr7CvACcCvwBeAPrLWD1toTwLeA3yptdx/wtLV2r7U2A3wFuNUYs3mOX8OCKB+fDxCqqyMfbaYpmAXg5Cs/wff9SpUnInJZlw18a+1PrLX/dnLaGJMAPgbsp9hVc6Rs9beA60qPt5Yvs9amga6y5UvK9AO4qfEcwYZlbG0YBiA31EeXVStfRBavKzpoa4xpAv4J+Dnwaml2umyVNBArPa6ftmz68iWl/ADuwHCGxniYX9q+ng42srNnPc+fXcXPXtCIHRFZvGY8Dt8YswX4HsVW+71AtLQoWrZaDBgrPU5NWzZ9+ZJTPj5/8qya4cYV1I8MMV4I8JNTAVa++Bq33bG9wpWKiLzbjFr4xpjbKLbqvwt81lqbsdYOAj0UD9pOeh/nunGOlC8zxsSAdZzfBbRkTZ5Vs7k5TqhxGRG3QNjN8/xLlkIhX+nyRETe5bItfGPMJuAZ4HettX86bfHfAl8zxhyk2IXzO8AjpWVPAHuNMXcALwMPA/utrY4rgfcm07Q0RQBoXrmKs2OD1DkFukcdvvW/XyDtxGhLxNixvX3qU4GISCXNpIX/ZaCB4lDMsbLbN4GvAm8Ah4FXKA7ffBTAWnsIeKA0fRa4Fvjc3L+EyigftRMOhwk2tDCUCzGWD9Hf3UNLYx0jqSyP7bJ0dCYrXK2IyAxa+Nbah4CHLrHKl0u3C227k+KbQNXZsb2dx3YVR+XEIkHCTcsZ7MmQCE0Q9nOMDpylcUUrUOz+UStfRCpNp1Z4j6aP2mlpjtHcEKEpWDxHfnqgF79QIBYJ0pucPlhJRGTh6WyZs1A+agfgT708J44MEyGHn/cY7u0m1NxKW2JJjkQVkSqjFv4cuvtDGyjEEmQKLr4PQ2cHGB3LsGN7e6VLExFR4M+lLesT/OvP3UwkHGA0HyLienw4MaD+exFZFNSlM8e2XrWC+z55A717vweAO9DNaN8tNLSurnBlIlLr1MKfB2bbteSb1gBQ8MHu+ZFOrCYiFafAnweu67DpozvwKZ4ff6z3FD1HdZEUEaksBf482bCxncDqa6amT7z0PIVctoIViUitU+DPo/ffcSd5NwzAxNgox3/xkwpXJCK1TIE/jxKJJhqu+fDUdO+hfWSGBipYkYjUMgX+PNt264fIRVsA8DwPu+dZHcAVkYpQ4M+zSDjImg99fGp66NQJ+o9VxRmiRWSJ0Tj8BRCoT/Cj1DWMjKZpCmY5+9yL/Nr6TQTCkUqXJiI1RC38edbRmeTxXZa6xEriwQLjhQAvdDfw42d+VOnSRKTGKPDn2eSVsRLNMcItK6eujPXioX5Gz3RWujwRqSEK/HnWm0wTixR7zpa3riAfilPnFBjKhTn6wg/wvVyFKxSRWqHAn2flV8ZyXYemVe1k/CBNwSwjySSdv3ixwhWKSK1Q4M+zHdvbSWU8UuM5fN/HDQbJRZZxbf0QAKcO/oJU3+kKVykitUCjdObZ5JWxdu/rojeZpi0R44ufvYWevd2QHsfzCuz+3rOcar6RvqGMLnwuIvNGgb8Apl8ZCyDu3kPnrsfpGw/xUk8dy0a7WLG2ferC5/fdbRT6IjKn1KVTIVdfvRZ304c4PNZM2M3jjybx0qPEoyHikSC793VVukQRqTIK/Aq65fYPk/QbqHMKFHyf4TMn8fOeLnwuIvNCgV9B8WiYdRvayfghAHLZHCPdJ0lnPF34XETmnAK/wj5z5xa8WMvUhc+Tg6MMnR3Qhc9FZM4p8Ctsy/oE//bzNxOKNTCaDxF189wc6GBd/USlSxORKqPAXwS2XrWc3/7iXfxK+wgfX95NIjBKx/PfpZDNVLo0EakiVzQs0xhzC/CMtba1NF0HjALl1+57yVp7V2n554FvAKuAPcD91tq+uSi82rQtr6ftQ3czuHcnAT9Hsj/JiRd/wMYdn8ZxnEqXJyJVYEaBb4xxgAeBb01bdB2QtNauvMA2W4HvAPcA+4BvAk8Cd86m4Gr2gW1X8aMzHyFwfA8ApzveoqH1ZVZs+0iFKxORajDTLp2vA18Cfn/a/JuAAxfZ5j7gaWvtXmttBvgKcKsxZvN7qrQGOI7DbXd+iFRiCwD5gs/xn71I+szxClcmItVgpoH/qLX2Joot9XI3Aq3GmIPGmF5jzD8YY9aUlm0Fpi7tZK1NA10UPxXIRUTqgtx09z2k6pYDkMnmsbufJjearHBlIrLUzSjwrbVnLrIoBfwU2AEYYBx4qrSsHpj+7aE0oAHml7GypZ51H/sk2UAUgOHhFMd+9H8p5DRyR0Teu1mdS8da+1D5tDHmIaDfGNNO8c0gOm2TGDA2m59ZK95vVjPQ/3EKB3+A6+fpPdNL5IV/YsOOz+C4gUqXJyJL0KyGZRpjfs8Yc03ZrHDpPkOxO8eUrRsD1lHWzSMX5zgOH/vI+8m0f3BqXtfRo/S+8jy+71ewMhFZqmY7Dn8b8EfGmGZjTDPwCPB9a20/8ATwq8aYO0rDNx8G9ltrO2b5M2tGIOByx10fZTixle5MlOf6VvLwM/088r+fp6NTffoicmVmG/gPAoPAMeAExfH4vwlgrT0EPAA8CpwFrgU+N8ufV3NikRBt132QF0fWMV4IEHdznDndzV9/91WFvohckSvqw7fWvgA0l00PAPdeYv2dwM73WpwU7e/oZ8XqVWT6TuIUMoTIkx/qYdeLR9jymx+tdHkiskTo1ApLQG8yzYpEjPpV6/Gc4pk13YLHO8dPkh242AAqEZHzKfCXgMkLoSeaY9S1riPvBJjwXWJk6HhuJ97I2UqXKCJLgAJ/CSi/EPqKlgayDavJFEJcWz/E2YEROp79P3j6YpaIXIYCfwmYvBB6YzzMwHCGTeuWc/MNG2mN5gDo6xvk2K6/xxsdrHClIrKY6SLmS8T0C6F7+QLPPR+g7tj/I+B79PQkcZ77ezbd9S8INiyrYKUisliphb9EBQMun9hxE6mNt1Nwit+87e4e4Nhzf6/uHRG5IAX+EhYKBrjr49sZWX/bVOj3dA9w7Id/R26ot8LVichio8Bf4iLhIL981y2MbbiNvFPsoevpG+LYs/+H7NlTFa5ORBYTBX4VKIb+B0ld9UvkS+P0e8+OYHf9I5lunUtfRIoU+FUiHApwzyduIr354+TcOgAGBlO8tespxo6/rhOuiYhG6VSTUDDAPTtu4PlwmMKbP6Iun2ZobIK3XtiF804/Px9aTt/gOG2JGDu2t5836kdEqp9a+FUmGHC567Zrqbvxn5MOFYdnHh0M8jd7euk/2UlLQ4iRVJbHdlmdfE2kxijwq5DrOtx+yyZaPvwphutWcXismaDjkRtNkjrzNrFggXgkyO59XZUuVUQWkAK/SjmOwwe3rWXLjl/lrN9MnVPA92F0JMVI13Ei/ji9yelXoBSRaqbAr3LXXNXC1ZvXMRxcjo8DQHp8gp7OkyxjFL+Qr3CFIrJQFPg14FMfu4qmlgSjkZV4BMgUXMay0J5+k94XnyKfHq10iSKyABT4NWDL+gT3/7NrMJtWkoqtJhgK8ZHmPhKBMY4e6eCdZx9j4syxSpcpIvNMwzJrxOTJ13zf59Cxft58aS/+6FvkCz6nziQZ+9F3WXftdTRvux03HKl0uSIyD9TCrzGO47Btcyu3f/Ju+tbcRjYQBWBodII3971G167HmOg9UdkiRWReKPBr1MqWOL/2yQ8TuOGfMxhZC0DOK/BOZy92106Sr+yikElVuEoRmUsK/BpWFwqw48ObufauT9G9/Ba80ikZkiMZ3nz1NU788G8Y73wD3y9UuFIRmQvqwxe2rFvG6hV38OIv1jJ69Bcsy5wi5xU4eTrJ4OAPWb3+EC3Xf4zw8rWVLlVEZkGBLwDUR0Pcc9v7OLpxJft+9jotA/upy6cZTWfpeOttlvd0s+rqLTRuvVVX1BJZohT4MsVxHLasW0Z720d5+fX19B7ex4r0cVw/T//gOIOvHaL1xNus3HoD8c03EYjEK12yiFwBBb68S7QuyJ23bKRnUyt7f3GMwJkDLBs/RVeqjmd76xl94ySrGk5y982r2fahmwlEGypdsojMgA7aykWtbInzmbu3Ye78JAfrb+XF4TWMFwLE3Rz9Y3n+6sfd7HnyCUYP7SGfHql0uSJyGWrhyyW5rsPWjS38eF8zjW0BMqMjBHNDRApZMgWPH3eGiGRfYvmRA6zY9D7qr7qO4LJVOI5T6dJFZJorCnxjzC3AM9ba1tJ0GPgz4LNAHvhja+3DZet/HvgGsArYA9xvre2bo9plAZ0dGmdtaz2F5fWcHV7GyOAg4YlhhnMBsrkCZ/rH6Eu+RsIeZvnq1TRefQN1q6/GCahNIbJYzOi/0RjjAA8C35q26OuAATYBTcCzxpjT1tq/McZsBb4D3APsA74JPAncOUe1ywJqS8QYSWWJR0O0JWIsb4pyqr+eQDrFWHiM+uwAXr5A3+A4/YPHaTzRRSLRSGLTViLthmDjCrX6RSpspn34Xwe+BPz+tPlfAP7AWjtorT1B8Q3ht0rL7gOettbutdZmgK8AtxpjNs++bFloO7a3k8p4pMZz+L5PJusRrQvzxc/ewurbP8OZNXeSjK6jgIsPDKeyvNN1lkM/3cvRZx6n98d/R/rtA/r2rkgFzfTz9qPW2q8aY+6YnGGMaabYVXOkbL23gOtKj7dSbNkDYK1NG2O6SsuPzqZoWXhb1ie4727D7n1d9CbTtCVifPr2TVPXxb3+6uUcP/0+3rBnyJ62JMZPEM6Pk80V6BtM0zf4DvETXTTV7yGxZh3xtVcTbttIIFpf4VcmUjtmFPjW2jMXmD35n1p+2aQ0ECtbPv2SSuXLZYmZPOPmhQQCLlvWLWPLumWcHbqKN471c+r4ceKjJ2ma6Mb186TGPVLjY3T3HyFuj9IQD9O8ai0NazcRbl1PoH6Zun1E5tFsjqhNfjaPls2LAWNly6Ocr3y5VKnlzVHu2L4O7wNrefv0MG+93cdo11Gax08Rzw7gA2PjOcbGc3SftUQ7jlMfC9HQ3Ezz2quoa11HaPka3FBdpV+KSFV5z4FvrR00xvRQPGh7ujT7fZzr4jlSWgaAMSYGrOP8LiCpYsGyVn9q/CqOnRri+Ds9jHe/TdNED/HsAA4+4xMe4xMe/YPjBDp7iEdfIR4N09C6ksaVawklVhNKrNJ5+kVmabZj5v4W+Jox5iDFLpzfAR4pLXsC2Fvq938ZeBjYb63tmOXPvKSOzuR5/cw7trdftBtCFk48GuL6zSuIhgP8cNDjtcFmAq7P2vAI6wNnqc/2E/A98gWfkVSWkVSW7rPHCNq3iUWCxCIh4onlNLauIryslWDTCoKNyzXsU+QKzPa/5avAHwGHKY74+UvgUQBr7SFjzAOl6TXAz4HPzfLnXVJHZ5LHdlnikSAtTRFGUlke22W5726j0F8Eyn8/G1Y3ks549KVDbL3+enzf50xXF8HRMzRM9BP1hgHw8oWpNwAGUjjHOomEg0TrgkTqQsSWLae+dRV1zcsJ1C8jUL8Mty6mYwGyoPKeh5fLkctm8bI5vFwWL5vFyxUfF3Ienpcjn8tS8HLkcx6FfI5CLkch71HwPPy8R9dQgYP9QYazQdrb2/jUnVvnNLsc3/fn7MnmkjFmA/DO7t27Wbt2Zqfl/Yudr0+NFZ+UGs/RGA/zpV+/fn4KlRm73O/H930GhjOc7hvjdE+SoTNdhNL9xHNJorlhHC78t+oA4VCAunDpFo0SaWqhvmU54cYEgVgjgTLBWDIAAApzSURBVFgjbrRBxwVqjO/7+F4xiCdvXrYUzKUwzudyxZtXvC94OfKeVwrjXCmMS/eFYjBTeszk9BzkaHcmyktDrYTdPHVOgWwwRqxt3RU3WE+dOsWOHTsANpaGy0+pqs/Dvck0LU3n9/PGIkF6k9MHC0klXO734zgOy5ujLG+Ocv2WFRQKW0iOZOg+m6K7f5jhnm68kX6iuWGiuSHq8sVxAz4wkcszkcsXhwoMjsOZJHCUYMAlHHQJhVzCwQChSIRwfSORhmaijU2EYvW4kThuXWzq5oTq9Alhnvm+D4U8fr4YmHnPI5edOBfG2RyeV2wl53M5vPJQ9kotY8+j4E22kHNTzzV5I18K6MXZpn2Xw2PNhN08Ebd4waF4PEosEmT3vq45a+VXVeCXfxt0Ujrj0ZbQSNDF4Ep/P6577g3guquXA5vITHj0DqbpHxzn7MAIY/09eMP91OVHqfPGqPPGCPje1HN4+QJevgATk3PGgUGgc+pnBAMOoYBLMOASCLgEg0HcugjBugjBuiihSPE+HIkSikQJ1EWKbwqBYNkthBMMTT0mEFwybxq+XyiGb6F0n8/jF/LFx36evJcnP9Ul4U2F77kuCu/8MJ4MYi831VVR8HLF8PXy+IVc6Wf4FPzibamE8iW5QQgEwQ1e4G8jiBsI4gRDuIEgbihEIBDEDYYIBEMEQkEmXhwkUR8kEAgSCIWIxmP4vj+nDdaqCvwd29t5bJcFii3HdMYjlfH49O2bKlyZwNz8fk72jEw7KH8tG9c0kRyZIDmcITkyzsjgEOmhAXIjSUK5McL5NOH8OKF8GpfzL9dYKPhkCz7Z3PTLOF589LDjFN8oAq6D6zg4joNbeuw6xeWO4+C4LrgBHNfFcUqPHRfHDeC7TvGxEwC3+Mbg4IAD/tQbxfR7Sp1aDviFYkj6fukSlD5+wS/N90u3AhQKU4/9gn+uZT0Z7oV8cTk+foGp+8LUc3CRjrSlo4BbPLg/GcCTgRwshrEbDOIGQnSnghzscxnOQKI+yE0bImxoixMIhQiWboFQuPg4HCIUChMMhwnVhQkEArju7E4+fFXXu7s857rBWlWBf7lvg0plzfb3c7mD8uf+MdYAxTAfTWenDvqOjE2QGh0lMzJMNjWMlx7F9TIECxOE8sX7YGHivE8IF+L7kM/75PNLPQory8el4AQoOC6+G8RxAxAIlVrCZWEcLLWGQ8XWsBsMTYVwYFoQh8JhQuHS47owoXCIUDBw2U9bHZ1JfrDLEm8Jkig1Rl7u9ti07eoFy4+FaLBWVeDDpb8NKpU3m9/P7n1dxCPBqRbQ5P3F+jhd16Gpvo6m+gsfqPV9n4lsnlQmRzpT/C7AeMYjPT5OtnTzMuPkJjIUshny2Qx+bgK3kGVgzOfYSJjRrEtjKMc19aOsio7j+nncgveuTxKLVXcmyuGxZoa98NTraIvm8B23FMguOAGcQPGTCW6g2C3hBnGDZV0UweC5MA6WwjkUJlgK51AoRDAcLraMS/fhULDYlRZ0OX5qiB9PNgSaF3449ZX+bc2HhWiwVl3gS/Wa64PyjuMQqQsSqQvS0jSzbXzf5813Bnh8l6VhTYDloSCpTI7DGY91N69nbVsD+byP5+UpFEr93/lCcTRHwSefL/aLky9wuj/Fwc4RhtMeTbEg718bY3WijqlOlMmO7VLj1PHBcfyp2p1i31Gpa8gF18F1XNxAaZ7r4jouTsAtdTm5BIIubiCAGwhysj/F/pe7iDWFWBcJksnmOTKR59qPb8asWzZ1TCPgzu+xiI7OJI9XeDj1YhnwMd8NVgW+LBmL4aC84zjs2X+axnh4qo5oJEhqPMdbnYP80vZ1M3qejs4kPzxoiTc2sqK1+PH9jT6PGz6wYcFCbte+12lpipa9jhDh8Rw/O9TNB7a0LkgNsDha14vhb2sh6BKHsmRMP0VzajxHKuOxY3v7FT1PR2eSv9j5Ov/1f73MX+x8nY7O5BVt35tME4uc31a60tZgecg5jkM8GiJeGoK3UObidVRLHXP1t7XYKfBlyZjs42yMhxkYztAYD1/xx/7JA78jqex53QdXEvptiRjpzPkHdq+0NbgYQm4uXke11DEXf1tLgbp0ZEmZbR/nXHQfzMVoisXQhbBYhjEvljpqYcCHWvhSU+aiZT0XrcHF0IWwWFq1i6WOWqAWvtSUuWpZz7Y1uFi+M7JYWrWLpY5qp8CXmrJYug9AIScLT106UlPUfSC1TC18qTlqWUutUgtfRKRGKPBFRGqEAl9EpEYo8EVEasRiPmgbAOjp6al0HSIiS0ZZZgamL1vMgb8K4N577610HSIiS9Eq4Hj5jMUc+K8AHwO6gXyFaxERWSoCFMP+lekLHL8qrh4sIiKXo4O2IiI1QoEvIlIjFPgiIjVCgS8iUiMU+CIiNUKBLyJSIxT4IiI1YjF/8eo9McZcDzwKbAPeBh6w1r7rCwhyYcaYB4BvAxNls78M/B3wZ8BnKX4R7o+ttQ8vfIVLgzHmFuAZa21raTrMJfafMebzwDcofmFmD3C/tbZvwQtfhC6wL+uAUSBbttpL1tq7Ssu1Ly+iqgK/9E/1PeBPgNuAXweeM8ast9aOVLS4peNG4I+stf+lfKYx5mHAAJuAJuBZY8xpa+3fVKDGRcsY4wAPAt+atujrXGT/GWO2At8B7gH2Ad8EngTuXLDCF6FL7MvrgKS1duUFttG+vIRq69K5AwhZa//EWpuz1j4JHAb+RWXLWlJuAg5cYP4XgD+w1g5aa09Q/Cf8rYUsbIn4OvAl4Penzb/U/rsPeNpau9damwG+AtxqjNm8QDUvVhfblxf7GwXty0uqtsDfCrw5bd5bFFsEchnGmADFrrDfNMacMcYcM8b8F2PMMoofj4+Ura79emGPWmtvoti6BMAY08yl99/W8mXW2jTQhfbvu/ZlyY1AqzHmoDGm1xjzD8aYNaVl2peXUG2BXw+kp81LA7EK1LIUraD4z/XXwEaK/c1fAv59aXn5vtV+vQBr7ZkLzK4v3V9s/+nv9gIusi8BUsBPgR0Uu8nGgadKy7QvL6Gq+vAp/iFEp82LAWMVqGXJsdb2ALeXzTpgjPlTiv2hcP6+1X6duVTp/mL7T3+3V8Ba+1D5tDHmIaDfGNOO9uUlVVsL/wjFd/xy7+P8j9JyEcaYa40xX582OwxkgB7O37farzNkrR3k0vvvvL9bY0wMWIf27wUZY37PGHNN2axw6T6D9uUlVVsL//8BjjHmtykOgft1in3ST11yK5k0BPwnY8wpiiMdPgD8B+DfUTz4/TVjzEGKH5t/B3ikUoUuQX/LxfffE8BeY8wdwMvAw8B+a21HJQpdArYB240xv1GafgT4vrW23xijfXkJVdXCt9ZmKXY//DqQBH4X+DVrbX9FC1sirLWngU9RHD0yAuwE/pu19h+BrwJvUAz+V0rLHq1QqUvRRfeftfYQ8EBp+ixwLfC5ypS5JDwIDALHgBMUx+P/JmhfXo4ugCIiUiOqqoUvIiIXp8AXEakRCnwRkRqhwBcRqREKfBGRGqHAFxGpEQp8EZEaocAXEakRCnwRkRrx/wHtfHMdYicbdgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plot(results.G, '-')\n",
    "plot(results2.G, '-')\n",
    "plot(data.glucose, 'bo')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The differences in `G` are less than 1%."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.970403647921524"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "diff = results.G - results2.G\n",
    "percent_diff = diff / results2.G * 100\n",
    "\n",
    "max(abs(percent_diff.dropna()))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Dropping pennies\n",
    "\n",
    "I'll start by getting the units we need from Pint."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "second"
      ],
      "text/latex": [
       "$\\mathrm{second}$"
      ],
      "text/plain": [
       "<Unit('second')>"
      ]
     },
     "execution_count": 19,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "m = UNITS.meter\n",
    "s = UNITS.second"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And defining the initial state."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>values</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>y</th>\n",
       "      <td>381 meter</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>v</th>\n",
       "      <td>0.0 meter / second</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "y             381 meter\n",
       "v    0.0 meter / second\n",
       "dtype: object"
      ]
     },
     "execution_count": 20,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "init = State(y=381 * m, \n",
    "             v=0 * m/s)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Acceleration due to gravity is about 9.8 m / s$^2$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "9.8 meter/second<sup>2</sup>"
      ],
      "text/latex": [
       "$9.8\\ \\frac{\\mathrm{meter}}{\\mathrm{second}^{2}}$"
      ],
      "text/plain": [
       "9.8 <Unit('meter / second ** 2')>"
      ]
     },
     "execution_count": 21,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "g = 9.8 * m/s**2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "When we call `odeint`, we need an array of timestamps where we want to compute the solution.\n",
    "\n",
    "I'll start with a duration of 10 seconds."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "10 second"
      ],
      "text/latex": [
       "$10\\ \\mathrm{second}$"
      ],
      "text/plain": [
       "10 <Unit('second')>"
      ]
     },
     "execution_count": 22,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "t_end = 10 * s"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we make a `System` object."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>values</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>init</th>\n",
       "      <td>y             381 meter\n",
       "v    0.0 meter / secon...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>g</th>\n",
       "      <td>9.8 meter / second ** 2</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>t_end</th>\n",
       "      <td>10 second</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "init     y             381 meter\n",
       "v    0.0 meter / secon...\n",
       "g                                  9.8 meter / second ** 2\n",
       "t_end                                            10 second\n",
       "dtype: object"
      ]
     },
     "execution_count": 23,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "system = System(init=init, g=g, t_end=t_end)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And define the slope function."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [],
   "source": [
    "def slope_func(state, t, system):\n",
    "    \"\"\"Compute derivatives of the state.\n",
    "    \n",
    "    state: position, velocity\n",
    "    t: time\n",
    "    system: System object containing `g`\n",
    "    \n",
    "    returns: derivatives of y and v\n",
    "    \"\"\"\n",
    "    y, v = state\n",
    "    g = system.g    \n",
    "\n",
    "    dydt = v\n",
    "    dvdt = -g\n",
    "    \n",
    "    return dydt, dvdt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "It's always a good idea to test the slope function with the initial conditions."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.0 meter / second\n",
      "-9.8 meter / second ** 2\n"
     ]
    }
   ],
   "source": [
    "dydt, dvdt = slope_func(system.init, 0, system)\n",
    "print(dydt)\n",
    "print(dvdt)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we're ready to call `run_euler`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "'Success'"
      ]
     },
     "execution_count": 26,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "system.set(dt=0.1*s)\n",
    "results, details = run_euler(system, slope_func, max_step=0.5)\n",
    "details.message"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>y</th>\n",
       "      <th>v</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0.0</th>\n",
       "      <td>381 meter</td>\n",
       "      <td>0.0 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>0.1</th>\n",
       "      <td>381.0 meter</td>\n",
       "      <td>-0.9800000000000001 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>0.2</th>\n",
       "      <td>380.902 meter</td>\n",
       "      <td>-1.9600000000000002 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>0.3</th>\n",
       "      <td>380.70599999999996 meter</td>\n",
       "      <td>-2.9400000000000004 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>0.4</th>\n",
       "      <td>380.412 meter</td>\n",
       "      <td>-3.9200000000000004 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>0.5</th>\n",
       "      <td>380.02 meter</td>\n",
       "      <td>-4.9 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>0.6</th>\n",
       "      <td>379.53 meter</td>\n",
       "      <td>-5.880000000000001 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>0.7</th>\n",
       "      <td>378.94199999999995 meter</td>\n",
       "      <td>-6.860000000000001 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>0.8</th>\n",
       "      <td>378.256 meter</td>\n",
       "      <td>-7.840000000000002 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>0.9</th>\n",
       "      <td>377.472 meter</td>\n",
       "      <td>-8.820000000000002 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1.0</th>\n",
       "      <td>376.59 meter</td>\n",
       "      <td>-9.800000000000002 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1.1</th>\n",
       "      <td>375.60999999999996 meter</td>\n",
       "      <td>-10.780000000000003 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1.2</th>\n",
       "      <td>374.532 meter</td>\n",
       "      <td>-11.760000000000003 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1.3</th>\n",
       "      <td>373.356 meter</td>\n",
       "      <td>-12.740000000000004 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1.4</th>\n",
       "      <td>372.082 meter</td>\n",
       "      <td>-13.720000000000004 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1.5</th>\n",
       "      <td>370.71 meter</td>\n",
       "      <td>-14.700000000000005 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1.6</th>\n",
       "      <td>369.23999999999995 meter</td>\n",
       "      <td>-15.680000000000005 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1.7</th>\n",
       "      <td>367.67199999999997 meter</td>\n",
       "      <td>-16.660000000000004 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1.8</th>\n",
       "      <td>366.006 meter</td>\n",
       "      <td>-17.640000000000004 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1.9</th>\n",
       "      <td>364.24199999999996 meter</td>\n",
       "      <td>-18.620000000000005 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2.0</th>\n",
       "      <td>362.37999999999994 meter</td>\n",
       "      <td>-19.600000000000005 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2.1</th>\n",
       "      <td>360.41999999999996 meter</td>\n",
       "      <td>-20.580000000000005 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2.2</th>\n",
       "      <td>358.36199999999997 meter</td>\n",
       "      <td>-21.560000000000006 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2.3</th>\n",
       "      <td>356.20599999999996 meter</td>\n",
       "      <td>-22.540000000000006 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2.4</th>\n",
       "      <td>353.95199999999994 meter</td>\n",
       "      <td>-23.520000000000007 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2.5</th>\n",
       "      <td>351.59999999999997 meter</td>\n",
       "      <td>-24.500000000000007 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2.6</th>\n",
       "      <td>349.15 meter</td>\n",
       "      <td>-25.480000000000008 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2.7</th>\n",
       "      <td>346.602 meter</td>\n",
       "      <td>-26.460000000000008 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2.8</th>\n",
       "      <td>343.95599999999996 meter</td>\n",
       "      <td>-27.44000000000001 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2.9</th>\n",
       "      <td>341.21199999999993 meter</td>\n",
       "      <td>-28.42000000000001 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>...</th>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7.1</th>\n",
       "      <td>137.4700000000001 meter</td>\n",
       "      <td>-69.57999999999993 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7.2</th>\n",
       "      <td>130.5120000000001 meter</td>\n",
       "      <td>-70.55999999999993 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7.3</th>\n",
       "      <td>123.45600000000012 meter</td>\n",
       "      <td>-71.53999999999994 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7.4</th>\n",
       "      <td>116.30200000000012 meter</td>\n",
       "      <td>-72.51999999999994 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7.5</th>\n",
       "      <td>109.05000000000013 meter</td>\n",
       "      <td>-73.49999999999994 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7.6</th>\n",
       "      <td>101.70000000000013 meter</td>\n",
       "      <td>-74.47999999999995 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7.7</th>\n",
       "      <td>94.25200000000014 meter</td>\n",
       "      <td>-75.45999999999995 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7.8</th>\n",
       "      <td>86.70600000000015 meter</td>\n",
       "      <td>-76.43999999999996 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7.9</th>\n",
       "      <td>79.06200000000015 meter</td>\n",
       "      <td>-77.41999999999996 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8.0</th>\n",
       "      <td>71.32000000000016 meter</td>\n",
       "      <td>-78.39999999999996 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8.1</th>\n",
       "      <td>63.48000000000017 meter</td>\n",
       "      <td>-79.37999999999997 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8.2</th>\n",
       "      <td>55.54200000000017 meter</td>\n",
       "      <td>-80.35999999999997 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8.3</th>\n",
       "      <td>47.50600000000017 meter</td>\n",
       "      <td>-81.33999999999997 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8.4</th>\n",
       "      <td>39.37200000000017 meter</td>\n",
       "      <td>-82.31999999999998 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8.5</th>\n",
       "      <td>31.14000000000017 meter</td>\n",
       "      <td>-83.29999999999998 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8.6</th>\n",
       "      <td>22.810000000000173 meter</td>\n",
       "      <td>-84.27999999999999 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8.7</th>\n",
       "      <td>14.382000000000174 meter</td>\n",
       "      <td>-85.25999999999999 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8.8</th>\n",
       "      <td>5.856000000000174 meter</td>\n",
       "      <td>-86.24 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8.9</th>\n",
       "      <td>-2.7679999999998266 meter</td>\n",
       "      <td>-87.22 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9.0</th>\n",
       "      <td>-11.489999999999826 meter</td>\n",
       "      <td>-88.2 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9.1</th>\n",
       "      <td>-20.309999999999825 meter</td>\n",
       "      <td>-89.18 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9.2</th>\n",
       "      <td>-29.227999999999824 meter</td>\n",
       "      <td>-90.16000000000001 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9.3</th>\n",
       "      <td>-38.24399999999983 meter</td>\n",
       "      <td>-91.14000000000001 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9.4</th>\n",
       "      <td>-47.357999999999834 meter</td>\n",
       "      <td>-92.12000000000002 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9.5</th>\n",
       "      <td>-56.56999999999984 meter</td>\n",
       "      <td>-93.10000000000002 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9.6</th>\n",
       "      <td>-65.87999999999984 meter</td>\n",
       "      <td>-94.08000000000003 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9.7</th>\n",
       "      <td>-75.28799999999984 meter</td>\n",
       "      <td>-95.06000000000003 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9.8</th>\n",
       "      <td>-84.79399999999984 meter</td>\n",
       "      <td>-96.04000000000003 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9.9</th>\n",
       "      <td>-94.39799999999984 meter</td>\n",
       "      <td>-97.02000000000004 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>10.0</th>\n",
       "      <td>-104.09999999999985 meter</td>\n",
       "      <td>-98.00000000000004 meter / second</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>101 rows × 2 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "                              y                                   v\n",
       "0.0                   381 meter                  0.0 meter / second\n",
       "0.1                 381.0 meter  -0.9800000000000001 meter / second\n",
       "0.2               380.902 meter  -1.9600000000000002 meter / second\n",
       "0.3    380.70599999999996 meter  -2.9400000000000004 meter / second\n",
       "0.4               380.412 meter  -3.9200000000000004 meter / second\n",
       "0.5                380.02 meter                 -4.9 meter / second\n",
       "0.6                379.53 meter   -5.880000000000001 meter / second\n",
       "0.7    378.94199999999995 meter   -6.860000000000001 meter / second\n",
       "0.8               378.256 meter   -7.840000000000002 meter / second\n",
       "0.9               377.472 meter   -8.820000000000002 meter / second\n",
       "1.0                376.59 meter   -9.800000000000002 meter / second\n",
       "1.1    375.60999999999996 meter  -10.780000000000003 meter / second\n",
       "1.2               374.532 meter  -11.760000000000003 meter / second\n",
       "1.3               373.356 meter  -12.740000000000004 meter / second\n",
       "1.4               372.082 meter  -13.720000000000004 meter / second\n",
       "1.5                370.71 meter  -14.700000000000005 meter / second\n",
       "1.6    369.23999999999995 meter  -15.680000000000005 meter / second\n",
       "1.7    367.67199999999997 meter  -16.660000000000004 meter / second\n",
       "1.8               366.006 meter  -17.640000000000004 meter / second\n",
       "1.9    364.24199999999996 meter  -18.620000000000005 meter / second\n",
       "2.0    362.37999999999994 meter  -19.600000000000005 meter / second\n",
       "2.1    360.41999999999996 meter  -20.580000000000005 meter / second\n",
       "2.2    358.36199999999997 meter  -21.560000000000006 meter / second\n",
       "2.3    356.20599999999996 meter  -22.540000000000006 meter / second\n",
       "2.4    353.95199999999994 meter  -23.520000000000007 meter / second\n",
       "2.5    351.59999999999997 meter  -24.500000000000007 meter / second\n",
       "2.6                349.15 meter  -25.480000000000008 meter / second\n",
       "2.7               346.602 meter  -26.460000000000008 meter / second\n",
       "2.8    343.95599999999996 meter   -27.44000000000001 meter / second\n",
       "2.9    341.21199999999993 meter   -28.42000000000001 meter / second\n",
       "...                         ...                                 ...\n",
       "7.1     137.4700000000001 meter   -69.57999999999993 meter / second\n",
       "7.2     130.5120000000001 meter   -70.55999999999993 meter / second\n",
       "7.3    123.45600000000012 meter   -71.53999999999994 meter / second\n",
       "7.4    116.30200000000012 meter   -72.51999999999994 meter / second\n",
       "7.5    109.05000000000013 meter   -73.49999999999994 meter / second\n",
       "7.6    101.70000000000013 meter   -74.47999999999995 meter / second\n",
       "7.7     94.25200000000014 meter   -75.45999999999995 meter / second\n",
       "7.8     86.70600000000015 meter   -76.43999999999996 meter / second\n",
       "7.9     79.06200000000015 meter   -77.41999999999996 meter / second\n",
       "8.0     71.32000000000016 meter   -78.39999999999996 meter / second\n",
       "8.1     63.48000000000017 meter   -79.37999999999997 meter / second\n",
       "8.2     55.54200000000017 meter   -80.35999999999997 meter / second\n",
       "8.3     47.50600000000017 meter   -81.33999999999997 meter / second\n",
       "8.4     39.37200000000017 meter   -82.31999999999998 meter / second\n",
       "8.5     31.14000000000017 meter   -83.29999999999998 meter / second\n",
       "8.6    22.810000000000173 meter   -84.27999999999999 meter / second\n",
       "8.7    14.382000000000174 meter   -85.25999999999999 meter / second\n",
       "8.8     5.856000000000174 meter               -86.24 meter / second\n",
       "8.9   -2.7679999999998266 meter               -87.22 meter / second\n",
       "9.0   -11.489999999999826 meter                -88.2 meter / second\n",
       "9.1   -20.309999999999825 meter               -89.18 meter / second\n",
       "9.2   -29.227999999999824 meter   -90.16000000000001 meter / second\n",
       "9.3    -38.24399999999983 meter   -91.14000000000001 meter / second\n",
       "9.4   -47.357999999999834 meter   -92.12000000000002 meter / second\n",
       "9.5    -56.56999999999984 meter   -93.10000000000002 meter / second\n",
       "9.6    -65.87999999999984 meter   -94.08000000000003 meter / second\n",
       "9.7    -75.28799999999984 meter   -95.06000000000003 meter / second\n",
       "9.8    -84.79399999999984 meter   -96.04000000000003 meter / second\n",
       "9.9    -94.39799999999984 meter   -97.02000000000004 meter / second\n",
       "10.0  -104.09999999999985 meter   -98.00000000000004 meter / second\n",
       "\n",
       "[101 rows x 2 columns]"
      ]
     },
     "execution_count": 27,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "results"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {},
   "outputs": [],
   "source": [
    "def crossings(series, value):\n",
    "    \"\"\"Find the labels where the series passes through value.\n",
    "\n",
    "    The labels in series must be increasing numerical values.\n",
    "\n",
    "    series: Series\n",
    "    value: number\n",
    "\n",
    "    returns: sequence of labels\n",
    "    \"\"\"\n",
    "    units = get_units(series.values[0])\n",
    "    values = magnitudes(series - value)\n",
    "    interp = InterpolatedUnivariateSpline(series.index, values)\n",
    "    return interp.roots()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([8.86802711])"
      ]
     },
     "execution_count": 29,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "t_crossings = crossings(results.y, 0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "'Success'"
      ]
     },
     "execution_count": 30,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "system.set(dt=0.1*s)\n",
    "results, details = run_ralston(system, slope_func, max_step=0.5)\n",
    "details.message"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([8.81788535])"
      ]
     },
     "execution_count": 31,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "t_crossings = crossings(results.y, 0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here are the results:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>y</th>\n",
       "      <th>v</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0.0</th>\n",
       "      <td>381 meter</td>\n",
       "      <td>0.0 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>0.1</th>\n",
       "      <td>380.951 meter</td>\n",
       "      <td>-0.9800000000000001 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>0.2</th>\n",
       "      <td>380.80400000000003 meter</td>\n",
       "      <td>-1.9600000000000002 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>0.3</th>\n",
       "      <td>380.559 meter</td>\n",
       "      <td>-2.9400000000000004 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>0.4</th>\n",
       "      <td>380.216 meter</td>\n",
       "      <td>-3.9200000000000004 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>0.5</th>\n",
       "      <td>379.77500000000003 meter</td>\n",
       "      <td>-4.9 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>0.6</th>\n",
       "      <td>379.23600000000005 meter</td>\n",
       "      <td>-5.880000000000001 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>0.7</th>\n",
       "      <td>378.59900000000005 meter</td>\n",
       "      <td>-6.860000000000001 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>0.8</th>\n",
       "      <td>377.86400000000003 meter</td>\n",
       "      <td>-7.840000000000002 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>0.9</th>\n",
       "      <td>377.031 meter</td>\n",
       "      <td>-8.820000000000002 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1.0</th>\n",
       "      <td>376.1 meter</td>\n",
       "      <td>-9.800000000000002 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1.1</th>\n",
       "      <td>375.071 meter</td>\n",
       "      <td>-10.780000000000003 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1.2</th>\n",
       "      <td>373.944 meter</td>\n",
       "      <td>-11.760000000000003 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1.3</th>\n",
       "      <td>372.719 meter</td>\n",
       "      <td>-12.740000000000004 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1.4</th>\n",
       "      <td>371.396 meter</td>\n",
       "      <td>-13.720000000000004 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1.5</th>\n",
       "      <td>369.975 meter</td>\n",
       "      <td>-14.700000000000005 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1.6</th>\n",
       "      <td>368.456 meter</td>\n",
       "      <td>-15.680000000000005 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1.7</th>\n",
       "      <td>366.839 meter</td>\n",
       "      <td>-16.660000000000004 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1.8</th>\n",
       "      <td>365.124 meter</td>\n",
       "      <td>-17.640000000000004 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1.9</th>\n",
       "      <td>363.31100000000004 meter</td>\n",
       "      <td>-18.620000000000005 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2.0</th>\n",
       "      <td>361.40000000000003 meter</td>\n",
       "      <td>-19.600000000000005 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2.1</th>\n",
       "      <td>359.391 meter</td>\n",
       "      <td>-20.580000000000005 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2.2</th>\n",
       "      <td>357.284 meter</td>\n",
       "      <td>-21.560000000000006 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2.3</th>\n",
       "      <td>355.079 meter</td>\n",
       "      <td>-22.540000000000006 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2.4</th>\n",
       "      <td>352.776 meter</td>\n",
       "      <td>-23.520000000000007 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2.5</th>\n",
       "      <td>350.375 meter</td>\n",
       "      <td>-24.500000000000007 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2.6</th>\n",
       "      <td>347.876 meter</td>\n",
       "      <td>-25.480000000000008 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2.7</th>\n",
       "      <td>345.279 meter</td>\n",
       "      <td>-26.460000000000008 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2.8</th>\n",
       "      <td>342.584 meter</td>\n",
       "      <td>-27.44000000000001 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2.9</th>\n",
       "      <td>339.791 meter</td>\n",
       "      <td>-28.42000000000001 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>...</th>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7.1</th>\n",
       "      <td>133.99100000000016 meter</td>\n",
       "      <td>-69.57999999999993 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7.2</th>\n",
       "      <td>126.98400000000017 meter</td>\n",
       "      <td>-70.55999999999993 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7.3</th>\n",
       "      <td>119.87900000000018 meter</td>\n",
       "      <td>-71.53999999999994 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7.4</th>\n",
       "      <td>112.67600000000019 meter</td>\n",
       "      <td>-72.51999999999994 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7.5</th>\n",
       "      <td>105.3750000000002 meter</td>\n",
       "      <td>-73.49999999999994 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7.6</th>\n",
       "      <td>97.9760000000002 meter</td>\n",
       "      <td>-74.47999999999995 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7.7</th>\n",
       "      <td>90.4790000000002 meter</td>\n",
       "      <td>-75.45999999999995 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7.8</th>\n",
       "      <td>82.8840000000002 meter</td>\n",
       "      <td>-76.43999999999996 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7.9</th>\n",
       "      <td>75.1910000000002 meter</td>\n",
       "      <td>-77.41999999999996 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8.0</th>\n",
       "      <td>67.4000000000002 meter</td>\n",
       "      <td>-78.39999999999996 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8.1</th>\n",
       "      <td>59.51100000000021 meter</td>\n",
       "      <td>-79.37999999999997 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8.2</th>\n",
       "      <td>51.524000000000214 meter</td>\n",
       "      <td>-80.35999999999997 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8.3</th>\n",
       "      <td>43.43900000000022 meter</td>\n",
       "      <td>-81.33999999999997 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8.4</th>\n",
       "      <td>35.25600000000022 meter</td>\n",
       "      <td>-82.31999999999998 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8.5</th>\n",
       "      <td>26.97500000000022 meter</td>\n",
       "      <td>-83.29999999999998 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8.6</th>\n",
       "      <td>18.596000000000224 meter</td>\n",
       "      <td>-84.27999999999999 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8.7</th>\n",
       "      <td>10.119000000000225 meter</td>\n",
       "      <td>-85.25999999999999 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8.8</th>\n",
       "      <td>1.5440000000002243 meter</td>\n",
       "      <td>-86.24 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8.9</th>\n",
       "      <td>-7.128999999999776 meter</td>\n",
       "      <td>-87.22 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9.0</th>\n",
       "      <td>-15.899999999999777 meter</td>\n",
       "      <td>-88.2 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9.1</th>\n",
       "      <td>-24.768999999999778 meter</td>\n",
       "      <td>-89.18 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9.2</th>\n",
       "      <td>-33.73599999999978 meter</td>\n",
       "      <td>-90.16000000000001 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9.3</th>\n",
       "      <td>-42.80099999999978 meter</td>\n",
       "      <td>-91.14000000000001 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9.4</th>\n",
       "      <td>-51.963999999999785 meter</td>\n",
       "      <td>-92.12000000000002 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9.5</th>\n",
       "      <td>-61.22499999999979 meter</td>\n",
       "      <td>-93.10000000000002 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9.6</th>\n",
       "      <td>-70.58399999999979 meter</td>\n",
       "      <td>-94.08000000000003 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9.7</th>\n",
       "      <td>-80.0409999999998 meter</td>\n",
       "      <td>-95.06000000000003 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9.8</th>\n",
       "      <td>-89.5959999999998 meter</td>\n",
       "      <td>-96.04000000000003 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9.9</th>\n",
       "      <td>-99.24899999999981 meter</td>\n",
       "      <td>-97.02000000000004 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>10.0</th>\n",
       "      <td>-108.99999999999982 meter</td>\n",
       "      <td>-98.00000000000004 meter / second</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>101 rows × 2 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "                              y                                   v\n",
       "0.0                   381 meter                  0.0 meter / second\n",
       "0.1               380.951 meter  -0.9800000000000001 meter / second\n",
       "0.2    380.80400000000003 meter  -1.9600000000000002 meter / second\n",
       "0.3               380.559 meter  -2.9400000000000004 meter / second\n",
       "0.4               380.216 meter  -3.9200000000000004 meter / second\n",
       "0.5    379.77500000000003 meter                 -4.9 meter / second\n",
       "0.6    379.23600000000005 meter   -5.880000000000001 meter / second\n",
       "0.7    378.59900000000005 meter   -6.860000000000001 meter / second\n",
       "0.8    377.86400000000003 meter   -7.840000000000002 meter / second\n",
       "0.9               377.031 meter   -8.820000000000002 meter / second\n",
       "1.0                 376.1 meter   -9.800000000000002 meter / second\n",
       "1.1               375.071 meter  -10.780000000000003 meter / second\n",
       "1.2               373.944 meter  -11.760000000000003 meter / second\n",
       "1.3               372.719 meter  -12.740000000000004 meter / second\n",
       "1.4               371.396 meter  -13.720000000000004 meter / second\n",
       "1.5               369.975 meter  -14.700000000000005 meter / second\n",
       "1.6               368.456 meter  -15.680000000000005 meter / second\n",
       "1.7               366.839 meter  -16.660000000000004 meter / second\n",
       "1.8               365.124 meter  -17.640000000000004 meter / second\n",
       "1.9    363.31100000000004 meter  -18.620000000000005 meter / second\n",
       "2.0    361.40000000000003 meter  -19.600000000000005 meter / second\n",
       "2.1               359.391 meter  -20.580000000000005 meter / second\n",
       "2.2               357.284 meter  -21.560000000000006 meter / second\n",
       "2.3               355.079 meter  -22.540000000000006 meter / second\n",
       "2.4               352.776 meter  -23.520000000000007 meter / second\n",
       "2.5               350.375 meter  -24.500000000000007 meter / second\n",
       "2.6               347.876 meter  -25.480000000000008 meter / second\n",
       "2.7               345.279 meter  -26.460000000000008 meter / second\n",
       "2.8               342.584 meter   -27.44000000000001 meter / second\n",
       "2.9               339.791 meter   -28.42000000000001 meter / second\n",
       "...                         ...                                 ...\n",
       "7.1    133.99100000000016 meter   -69.57999999999993 meter / second\n",
       "7.2    126.98400000000017 meter   -70.55999999999993 meter / second\n",
       "7.3    119.87900000000018 meter   -71.53999999999994 meter / second\n",
       "7.4    112.67600000000019 meter   -72.51999999999994 meter / second\n",
       "7.5     105.3750000000002 meter   -73.49999999999994 meter / second\n",
       "7.6      97.9760000000002 meter   -74.47999999999995 meter / second\n",
       "7.7      90.4790000000002 meter   -75.45999999999995 meter / second\n",
       "7.8      82.8840000000002 meter   -76.43999999999996 meter / second\n",
       "7.9      75.1910000000002 meter   -77.41999999999996 meter / second\n",
       "8.0      67.4000000000002 meter   -78.39999999999996 meter / second\n",
       "8.1     59.51100000000021 meter   -79.37999999999997 meter / second\n",
       "8.2    51.524000000000214 meter   -80.35999999999997 meter / second\n",
       "8.3     43.43900000000022 meter   -81.33999999999997 meter / second\n",
       "8.4     35.25600000000022 meter   -82.31999999999998 meter / second\n",
       "8.5     26.97500000000022 meter   -83.29999999999998 meter / second\n",
       "8.6    18.596000000000224 meter   -84.27999999999999 meter / second\n",
       "8.7    10.119000000000225 meter   -85.25999999999999 meter / second\n",
       "8.8    1.5440000000002243 meter               -86.24 meter / second\n",
       "8.9    -7.128999999999776 meter               -87.22 meter / second\n",
       "9.0   -15.899999999999777 meter                -88.2 meter / second\n",
       "9.1   -24.768999999999778 meter               -89.18 meter / second\n",
       "9.2    -33.73599999999978 meter   -90.16000000000001 meter / second\n",
       "9.3    -42.80099999999978 meter   -91.14000000000001 meter / second\n",
       "9.4   -51.963999999999785 meter   -92.12000000000002 meter / second\n",
       "9.5    -61.22499999999979 meter   -93.10000000000002 meter / second\n",
       "9.6    -70.58399999999979 meter   -94.08000000000003 meter / second\n",
       "9.7     -80.0409999999998 meter   -95.06000000000003 meter / second\n",
       "9.8     -89.5959999999998 meter   -96.04000000000003 meter / second\n",
       "9.9    -99.24899999999981 meter   -97.02000000000004 meter / second\n",
       "10.0  -108.99999999999982 meter   -98.00000000000004 meter / second\n",
       "\n",
       "[101 rows x 2 columns]"
      ]
     },
     "execution_count": 32,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "results"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And here's position as a function of time:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Saving figure to file figs/chap09-fig01.pdf\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAagAAAEYCAYAAAAJeGK1AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjAsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+17YcXAAAgAElEQVR4nO3deXTc9X3v/+fMaLSMVlu7LdmSvHy8YAPewICxwWEPGAhZoZBf0obm5P5uml/p/ZXmhDS9paQloUmzXEJvTvtrWkKbS0hYAhgMNjY2YPC+6ONNsi3ZWqzF2teZ3x/f0XgkvAhbs2j0epwzR/PdZt46tuY9n8/383l/XIFAABERkXjjjnUAIiIiZ6MEJSIicUkJSkRE4pISlIiIxCUlKBERiUtKUCIiEpeSYh3ASMaYHGAX8Ji19l+NMcnAT4H7gEHgKWvtE2Hnfw74O6AY2AB82VrbEP3IRURkLMVdggKeBqaGbX8PMMAMIBt4zRhTa639N2PMPOCXwG3Ah8DfA88BN472zYwxKcBS4CROAhQRkejx4DQwtlpre8MPxFWCMsY8BGQBu8N2P4TTKmoBWowxPwAeBv4NeAB4yVq7KXj9o8FzZllrD47ybZcCG8fqdxARkYuyAtgUviNuEpQxphz4LnAN8FpwXw5OZt0XdmolsCD4fB5OywkAa22XMeZ48PhoE9RJgP/4j/+gqKjoUn4FERH5hOrq6rj//vsh+FkcLi4SlDHGA/w78Ii1ts4YM3QoI/izK+z0LsAXdjz82MjjozEIUFRURElJyScJW0RExs7HbrHEyyi+7wDWWvvbEfs7gz/Twvb5gI6w42kMF35cRETGqbhoQQFfAKYYY+4NbmcCPweWAXU4gyRqg8fmcKbLb1/wGADGGB8wjeFdgiIiMg7FRYKy1s4J3zbG7AB+FBxm3gF81xizC6dL7xHgx8FTnwU2GWNWAVuAJ4Dt1toDUQteREQiIi4S1AU8BvwQ2IvTJfkMzlB0rLW7jTFf4czQ9PeBz8YoThGRS+b3+6mpqaGzs/PCJ48DXq+XgoICsrKyPvG1cZmgrLVXhD3vAb4RfJzt3OeB56MUmohIRJ06dQqXy4UxBrc7XoYJXJxAIEB3dze1tc4dmk+apOIyQY0XbZ197Dl8ip6+QVwucAEulwu323l43C6SPG7cbhdJHhdejwdvkhtvkptkr4dkr5sUr4eUZA9JHjculyvWv5KIxFhraytlZWXjPjmB83no8/mYOnUqJ06cUIKKpvXbjnOsrn1MXsvtdpGanERqsoe0lKTQIz3Niy81ifRUL+lpXjLSvKQke5TMRBLU4OAgXq831mGMqbS0NPr7+z/xdUpQl8CXMnb/ifz+AF09/XT1XPgfMcnjJsPnJdOXTKYvmax055GdkUJ2ejKpKfpnFRnPEu0L6MX+PvokuwQ3LC6hYmo23b0DBAIBAjh9rn5/AL8fBv1+BgYDzs8BPwODfvoHnEffgJ++/kH6+gfp6RtkYNA/6vcdGPTT2t5La3vvWY+nJieRnZHMpMwUcjJTmZyVyqSsFLLTU3C7E+s/vogkLiWoS+DxuKmYmj0mrzUw6Kenb5Ce3gG6hx49A3QGW1WdPQN0dvfT0d1PX//5a9r29A3Q0zxAffPwIhset4vJWU7Cys1OIzcnlfycNHypidWdICKJQQkqTiR53GSkuclIu3Cy6O0fpKOrj7bOPtqDP9s6+zjd0cfpjt5ztsYG/QEaW7tpbO0GWkL7faleCialkZ+TRsFkH4WTfUpaIhJzSlDjUIrXQ0p2GrnZI6s8OV2Mnd39tLT30trRS2tbL83tPbS09dDRffb7W109/VSf7Kf6ZFtoX0aal8LcdAon+yjOTSd/UhpJnvE/qkhExg8lqATjcrnI8CWT4UumtDBz2LGevgGa23poOt1DU2s3p4I/+8/S4uro7qejppXDNa2AM8qwYJKP4rx0puSlU5ybrsEYIlGw3Tbwwb46+gdGf5/6UniT3CybV8SVpuC8533nO9+hp6eHJ598MrTvpptu4tFHH+XGG0e9JN956RNmAklNTmJKXgZT8jJC+/z+AK0dvTS2dNHQ0h36ObKb0O8PUNfUSV1TJ9utkwhzs1OZmpfB1IIMpuSnk5qs/04iY23HgcaoJSeA/gE/Ow40XjBBrVmzhq997Wv09vaSkpLCzp07aW9vZ8WKFWMWiz5RJjh32MAJM93Z5/cHaG7rob65i7qmTk42dX5sxGAgEOBUazenWrvZeagRl8tFfk4aJQUZlBZmMiUvHY+6BEUu2RWz86Pegrpidv4Fz1u8eDE5OTmsX7+eW265hZdffpnbb799TOdwKUHJx7jdLvJy0sjLSWN+RS7g3Keqb+7iRGMnJ0510NjSjT8QCF0TCARoaOmioaWLbbYBr8fNlPwMphdnMr0oi+yMlFj9OiLj2pWm4IKtmVhwuVzceeedvPLKK9x00028+uqr/OxnPxvT91CCklHxpXopn5JN+RRnWH1f/yAnmzqpbeigtrGDhpZuAmEJq3/Qz9G6No7WtQG15GSkUDYli+lFWUzJz8Cj+Vgi496aNWu499572bx5Mz6fj8svv3xMX18JSi5KstfD9CIn4YAzAONEYyfH69s5Xt9Oa8fwLsHWjl52HGhkx4HG0LXlU7KYXpxFitcTi19BRC5RRUUFM2bM4Mknn+Suu+4a89dXgpIxkZqcRMXU7NDE5dMdvRyra+dYXRs1DR3DRgr29Q9y8HgLB4+34Ha7KCnIYMbUHMqnZGn+lcg4s2bNGh5//HF+8pOfjPlrK0FJRGRnpLBgZgoLZuYxMOintrGDoyfbqD7ZRltnX+g8vz8QTGTtrN/mYmp+OjNLcqiYmq1kJTIOFBcXs2jRIqZNmzbmr60EJRGX5HGHugNXXBGg6XQPVSdOc+TEaRpbukPnBQIBaho6qGnoYMP2WkoKMphV6iQrDWEXiS/t7e3U1tbyzDPP8MUvfjEi76G/eokql+vMCMGl84po6+yjqvY0h2tbOdnUFRpoEQgEQvezNmyrYXpxFrOnTaKsOEsVLUTiQFVVFQ8++CDXX389a9asich7KEFJTGWlJ3P57Hwun51PZ3c/h2tbOXS8lROnzix3PegPcKT2NEdqT5Pi9TCzNAczbRLFeekJtyyByHixcOFCduzYEdH3UIKSuJGe5mXhzHwWzsyno6uPQzWtHDjWSkPLmarsvf2D7D3SxN4jTWRnpDC3bDJm+iQyfckxjFxkbAUCgYT68hU+BeWTUIKSuJThS+aK2QVcMbuAlvYeDh5rpfJo87ABFqc7enlvz0ne31tHaUEG88pzKZ+SpQoWMq55PB76+/tJTk6cL13d3d0XVWFCCUri3qTMVJbNL2LpvELqmrqoPNrMoeOt9AbXxQoEAhyrb+dYfTtpKUnMmT6ZeRWTmZSZGuPIRT65nJwc6uvrmTp1Km73+P6yFQgE6O7upra2lsLCwk98vRKUjBsul4vivHSK89JZccVUjtSeZn91MzUNHaEuhO7eAbYfaGD7gQam5mdw2YxcKqZkq1Ul40ZeXh41NTVYa2Mdypjwer0UFhaSlZX1ia9VgpJxKcnjZva0ScyeNom2zj4qq5vZX91Me9eZLsDaRqcMky/Vy/zyycyvyCVD96okzrnd7ojMKRqPlKBk3MtKT2bZ/CKWzC3keH07e440UX2yLdSq6urpZ+v+ej6qbKBiajYLZ+ZpBKDIOKAEJQnD7XYxvdip79fR1ce+qmb2Hmmis8dZSdgfCHCoppVDNa3k56Rx+ax8ZpXmqPtPJE4pQUlCyvA5rarFcwupOnGa3YdOUdvYETre2NrNm1uPsXn3SS6bkcuCGXmkaYVgkbiiv0hJaB63i5klOcwsyeFUaze7D5/CHm0JrRjc1dPPB3vr2FbZwJzpk7h8dr5G/4nECSUomTDyctK4YXEpyy8rZm9VE7sPnaKj2+n+Gxj0s+dIE3urmqmYksWiOYUUTvbFOGKRiU0JSiac1JQkFs8p5IrZBRyuaWX7gYZQ0dpAIMDh2tMcrj3N1PwMFs0pYFphpgZUiMSAEpRMWB63i9nTJjGrNIfaxg6228bgCsCOoWHq+ZPSWDynkBlTs5WoRKJICUomPJfLRUlBJiUFmTSd7ma7beDAsVb8wWHqjS3dvLalmslZqSyZW8jMkhzcWrJeJOI0vlYkTG52Gp9aNp0HbpvLwpl5w5b2aG7rYe37R3n29Uoqq5vx+y+uAKaIjE5ctaCMMZ8G/g4oBxqAf7DW/sIYkwz8FLgPGASestY+EXbd54LXFQMbgC9baxuiHb8kjqz0ZK6/soQlcwvZebCR3Yeb6AvW/mvt6OXNrcf4cH89S+YVMrt0klpUIhEQNy0oY0wx8H+A/9damwl8FviRMWYR8D3AADOApcBDxpgHg9fNA34JfBnIBQ4Cz0X9F5CE5Ev1snzBFB68bS7L5hWRkuwJHWvt6OXND47x7OuVHDjWohaVyBiLmwRlrT0J5FtrXzXGuHGSzQDQDjwEPG6tbbHWVgM/AB4OXvoA8JK1dpO1tgd4FLjWGDMr6r+EJKzUlCSWzS/iwdvncfVlxR9LVGvfP8p/vmE5XNN60WvfiMhwcZOgAKy17cYYH9ALrAV+BjTidN3tCzu1ElgQfD4v/Ji1tgs4HnZcZMykeD0smVvIg7fP46r5w1tUTW09vLqlmt+sO8jx+nYlKpFLFFf3oIJ6gHRgIfAHoDu4vyvsnC5gaBZlxohjI4+LjLkUr4el84pYMDOPnQca2XGwkf4BpzpFQ0sXv3/nMCUFGSxfMEUTfkUuUtwlKGutH+gDPjTGPAMsCR5KCzvNBwwVVusccWzkcZGISU1O4qrLilk4K59ttoHdh06FyijVNHTwm3UHmFGSw/LLisnJTIlxtCLjS9x08RljVhpjPhqxOwVoAepwBkkMmcOZbr194ceCXYTTGN4lKBJRaSlJXLtwCg/cNpf5Fbm4wyb0Hq5p5dnXK9mwrYauYGV1EbmweGpB7QCmGmP+H+DHwFXAV4F7cBLUd40xu3C69B4JngPwLLDJGLMK2AI8AWy31h6IbvgikJHm5YbFpVw5u4D39pzkUE0r4Cz1sfvwKeyxFhaZAq6YnT9sjpWIfFzc/IVYa08DtwP3As3AM8AfW2s3AI8Be4C9wFbgeeDp4HW7ga8Et08B83GGqIvETE5mCrcuL+Ozq2czNT8jtL+vf5D39pzk31/djz3arIEUIufhmuh/IMaYMqBq3bp1lJSUxDocSUCBQICjde1s3nWC5raeYccKJ/u47vKpFOelxyg6kdiqqalh9erVAOXBaUQhcdOCEklULpeLsuIsvnCTYdWikmELI9Y3d/H82wdZ+/5ROrr6YhilSPyJp3tQIgnN7XZx2Yw8Zk+bxEeV9ew40MhgsPrEgWMtVNWe5so5BSwyBbo/JYJaUCJRl+z1sHzBFO6/dS4zS3JC+/sH/Xywt45nX69URQoRlKBEYiYrPZlbl5dxz6qZ5OWcmcrX1tnHq1uqeWnjEVrae879AiIJTglKJMam5mfwudWzWbWohNTkM73ux+rb+fVay5bdJ0NVKkQmEiUokTgwdH/qgVvnsGBGXmjlXr8/wEeV9Tz7eiVVJ07HOEqR6FKCEokjqSlJrFxUwudWz6Y498zQ8/auPl55t4pXNh2hrVOj/WRiUIISiUP5k9K494aZrF4ybdiw9KqTbfz69Uq2VTaERgCKJColKJE45XK5mFs+mftvmcNlFbmh/f2DfjbvPsFv1h2grqkzhhGKRJYSlEicS01JYtXiUu67cdaw0X6nWrt5/u1DvLO9JrQcvUgiUYISGSeKctP53OrZXLNwCt7gRN5AIMCuQ6d49vVKqk+2xThCkbGlBCUyjrjdLhaZAr5ws2FaUWZof0d3Py9vOsLr7x3Vkh6SMJSgRMah7IwU7ryugpuvmj5sEMXB4y38eq3lwLEWVaKQcU8JSmSccrlczJ42iS/dMoc50yeF9nf3DrD2/aP8YXM1Hd1qTcn4pQQlMs6lpSTxqWXTuXNFBRlp3tD+qhOn+fXaSiq17pSMU0pQIglielEWX7plDpfNyAvt6+0b5M0PjvGHd6vUmpJxRwlKJIEkez2sWlTCPatmkpWeHNpfdbKNX6+t1L0pGVeUoEQS0NT8DL54s2HBiNbU2veP8ppG+sk4oQQlkqC8SR5WLirh7pUzhrWmDte08uu1VsVnJe4pQYkkuJKCTL5wk2F+WLmk7t4BXnm3irc+PK4qFBK3lKBEJoBkr4cbFpdy53UVpKeeGem3r6qJ596wnDylmn4Sf5SgRCaQ6cVZfPFmM2yp+bbOPn67/hDv7zmpCukSV5SgRCaY1JQkbrl6Ojctm0aK1wM4Nf227q/nt28fpLW9N8YRijiUoEQmIJfLhZk+mS/cbJianxHaX9/cxX++aams1uReiT0lKJEJLNOXzJrrZ3DNgim43c4y8/0Dft7ceoy17x+jp28gxhHKRJZ04VMcxphCYDFQAAwCdcA2a21ThGITkShwu10smlNASWEGa98/GuriO3i8hfrmTm6+ajpFYcvPi0TLeROUMSYJ+BLwZ8DlQB/QAniAycFz3gd+DjxnrfVHNFoRiZiCST4+/6nZbNxxgn1VzvfOts4+fvv2Ia66rIhFpgCXyxXjKGUiOWcXnzFmJbALeBD4JTAb8Flrp1hrC4Fk4ErgWeC/AZXGmFURj1hEIsab5OHGJaXcurwsNIDCHwiwZfdJXtx4RBUoJKrO14L6c+Dz1trdZztorQ0Ae4KPnxtjrgT+Blg/1kGKSHTNLMmhYJKPte8fpa7JmSN1vL6d5944wE3LplFamHmBVxC5dOdMUNbauz7JC1lrtwN3XnJEIhIXstKTuWfVTD7YW8c220AgEKCrp58XNx5h6dxClswtDA2sEImETzJIwgeUAykjj1lrt41lUCISHzxuF8sXFFNS4Ayg6O4dIBAI8MG+Ok6c6uDmq6bjC6tMITKWRjXM3BjzANCAc0/qwxGPrRGLTkTiQmmhU8+vpODMnKmahg6ee+MAtY0dMYxMEtloW1BP4AyUeAroiUQgxpibgO8Ds3CS4ZPW2l8YY5KBnwL34Qxvf8pa+0TYdZ8D/g4oBjYAX7bWNkQiRpGJLD3Ny10rZvDh/nq27q8Pdfn9bsNhrppfxOI5GuUnY2u0CSoL+Km19mgkgjDGlALPAw8Bv8eZb/W6MaYaWAUYYAaQDbxmjKm11v6bMWYeTuK8Dac19/fAc8CNkYhTZKJzu10sm19EcV76sC6/9/acpL6pk9XLppGaPOo7ByLnNdpKEr8CvhzBOMqAZ621L1hr/dbarTijAa/FSVqPW2tbrLXVwA+Ah4PXPQC8ZK3dZK3tAR4FrjXGzIpgrCITXmlhJp+/yTAl70yXX9XJNv7rzQM0NHfFMDJJJKP9qvMksM0Ycz9QDQybkGutvaQWi7V2I7BxaNsYMxlYgZMYi4F9YadXAguCz+fhtJyGXqfLGHM8ePzgpcQkIueXkeZlzcoZvLfnJNuDveptnX08//ZBVi4qYV557gVeQeT8RpugfgV0AK8AEf16ZIzJBl4E3gc+Cu4Of88uwBd8nnGWeMKPi0gEedwurl04haLJPtYFFz8c9Ad468Pj1Dd3cf0VU/F4VPJTLs5oE9RS4Cpr7a5IBmOMmY1zD2ofcD+QFjyUFnaaDydZAnSOODbyuIhEwYySHHKz03h1cxVNbc44qr1HmjjV2s1ty8vI8CVf4BVEPm60X20skHPBsy6BMeZ6nFbT74D7rLU91toWnKK0JuzUOZzp8tsXfiw4V2saw7sERSQKcjJTuG/1LGaVTgrtc5bv0FB0uTifZJj5vxpjfgocBoYV5LLW/uFSgjDGzABeBr5trf3JiMO/Ar5rjNmF06X3CPDj4LFngU3BGoBbgnFut9YeuJR4ROTieJM83HzVNIpyfby78wT+QIDu3gF+v+Ew1y6cwsJZeRqKLqM22gT16+DPH5zlWACnuvml+AaQCTxhjHkibP/PgMeAHwJ7cVp8zwBPA1hrdxtjvhLcnorTAvvsJcYiIpfA5XJx+ax88nLSeG1LNd29A/gDATburKWxtZtVi0tI0n0pGQXXRF810xhTBlStW7eOkpKSWIcjklA6uvp4dUs19WFDzwsn+7jtmnIy0lQiSaCmpobVq1cDlAenEoVcaLmNT8QYowmyIhKS4XMKzs4tmxzaV9/cxX+9eSBUJV3kXM7XxfctY8xfAv8EvGmtPetCMMFFDT+NsyZUF/DWmEcpIuNWksfNjUtKyZ+UxqYdzn2prp5+Xlh/iFWLSplbPvnCLyIT0vmW27jbGHMPTn286caY9Tj3gU4BLiAfZ5Xd5cAx4H9aa/9PxCMWkXHH5XKxcGY+k7NSeW3LUXr6Bhj0B1j34TGa2rq5ZsEULd0hH3PeQRLW2heAF4Kj5G7HSUaFOJUk6nAm0j4RrAQhInJeJQWZfHb1LP7w7pn5UjsONNLc1sMtV59ZxVcERjmKz1q7Hq2UKyJjIDsjhc/cOIs3tx7jSO1pAI7VtfP8Wwe549pysjM+tuScTFAa6ykiUZfs9XDb8jKWzC0M7Wtu6+E36w5qUq+EKEGJSEy4XC6uvqyYm6+ajid4/6mnz5nUu6+qKcbRSTxQghKRmJo9bRL3rJoZWjreH3CKzW7edYKJPk9zolOCEpGYK8pN53OrZ5GXc6b28zbbwGtbqukf8J/7Qkloo1760hhTACwEvDjDzEMutRafiEiGL5l7V83kjfePUnWyDYDDtafp2HCIO64tD7WwZOIYVYIyxnwV+DlOchppLGrxiYg4gyeuKWfz7hPsONAIOJUnfrPuIJ++rpzc7JGr60giG20L6i+AfwYetda2RzAeEZng3G4X110+leyMFN7ZXksgEKC9q4/n3z7EbcvLKC3MjHWIEiWjvQdVCvxYyUlEomXBjDw+fW053iTnY6qvf5CXNh5hf1VzjCOTaBltgloLrI5kICIiI00vzuIzN8wKVT73B5zySB/srdMIvwlgtF18O4GnjDF3AQeAvvCD1tr/MdaBiYgA5OWkcd+Ns3j53SpOtXYD8MG+Otq7+li1uDQ0h0oSz2gT1EqcxQDTcArEhtPXGBGJqKERfq9tqeZYvXOnYX91M53d/dy6vIxk1fBLSKOtxXdDpAMRETmfZK+HO66rYP1Hx9lf7dyHOlbfzgvrD/Hp6ypI1wKICeeTzIMqxFnzaT7Ovav9wD9ba49EKDYRkWE8bhc3Likl05fMB/vqAGhs7eb5tw9y54oKJmWmxjhCGUujGiRhjFmGc+/pHpz1oBpxFincZYxZErnwRESGc7lcLJtfxI1LSnG7nPtPbZ19PP/WIa3Sm2BG24L6IfBr4OvW2tA9J2PMT4EnAXUBikhUzSvPJT3V65RDGvTT0zfA7zYc5tblZZQVZ8U6PBkDox1mvgT4x/DkFPQTYOnYhiQiMjrTi7NYs3IGqcnOd+2BQT9/eLeKymrNlUoEo01QJ4Gys+yvADR5V0Ripig3nc/cOJOs9GTAmSv15tZjbLMNMY5MLtVou/h+BTxjjPkz4L3gvuXAPwaPiYjEzKTMVD5zwyxe2nQkNFdq864TdPcMcM3CYlwuzZUaj0bbgnocp5rEfwE1QC3OPanfAN+OTGgiIqOXnubl7pUzmJKXEdq3/UAD67Yex+/XdM3xaLTzoPqAPzHGPAIYoBs4ZK3tjmRwIiKfRGpyEnddX8Eb7x/lcO1pACqPNtPbP8gtV08nyaMl8MaTcyYoY8ztwBvW2v7g85FKjTGA1oMSkfiR5HFzy9VlrN92nH3BwrJVJ07z0sYj3HFtuapOjCPna0G9DBQBDcHn56L1oEQkrrjdLm5YXEpqclJosERtYwe/23CYO1dUkJYy6hoFEkPn/Fey1rrP9lxEZDxwuVxcs3AKqSlJbN51AoCGli5eWH+Iu1ZUkOFLjnGEciGjrSTxljEm5yz7840xH419WCIiY2ORKeCGxaWhkXzNbT38dv0hWtt7YxyZXMj57kGtAuYFN1cCDxtjRs55mgvMiExoIiJjY35FLsleN298cAy/P0BbZx+/XX+INddXaBn5OHa+jtgm4BHAFXx8AxgMOx4AOoA/j1h0IiJjZFbpJJK9Hl7dXM3AoJ+unn5eWH+Yu1ZUUDDZF+vw5CzOdw9qN06lCIwxbwP3WmtbIh1QsDDty9baguB2MvBT4D6cBPmUtfaJsPM/B/wdUAxsAL5sraaQi8jHTS/K4q4VFbz8bhV9/YNO/b53DvPpa8uZkp9x4ReQqDrnPShjTPhXijuAXmOM72yPsQjEGOMyxvwxzoTg8LuX38OZezUDp+7fQ8aYB4PXzAN+CXwZyAUOAs+NRTwikpim5Gew5voz9fv6+gd5ceMRjteralu8Od8giXZjTEHweQdOzb2Rj6H9Y+F7wNeBvx2x/yHgcWtti7W2GvgB8HDw2APAS9baTdbaHuBR4FpjzKwxiklEElDhZB/3rJqBL9VZ5HBg0M/Lm45QdeJ0jCOTcOdLUDcCQyWBbwhuj3wM7R8LT1trFwMfDu0IjhwsBvaFnVcJLAg+nxd+zFrbBRwPOy4icla52Wncs2oGGcGVeAf9AV7dXM2hmtYYRyZDzncPasPZnkPovtBC4IC1tm0sArHWnjjL7qFO4a6wfV2AL+x4F8OFHxcROadJmance8MsfrfhEG2dffgDAV5/7yj+ZQFmT5sU6/AmvNHOg5ppjNlgjLk6eM/pg+DjqDHm6gjGN7Q8Zvg4UB9O1+LQ8ZFjRMOPi4icV1Z6MvfeMIuczBQAAoEAb3xwjP1VWlMq1kZbIeInOPeaqoE/AkpwBi78L+CpiEQGBEcN1gXfa8gcznTr7Qs/Fkye0xjeJSgicl4ZaV7uXTWT3KxUwElS6z48xt4jTTGObGIbbYJaAXzLWlsH3A28Yq09CPwzcEWkggv6FfBdY0yeMaYMZ27W0BpUzwJrjDGrjDEpwBPAdmvtgQjHJCIJxpfq5e5VM8nPOdMp8/ZHx9l96FQMo5rYRpugegCvMSYdp6rEq8H9RUCkh708BuwB9gJbgeeBpyE0V+srwe1TwHzgsxGOR0QSVFpKEpUDlH8AABN0SURBVGuun0HBpDO3sTdsr2HXocYYRjVxjbak7+s4raV2nEEILxljVgM/Bl4cy4CsteuBnLDtHpwqFt84x/nP4yQtEZFLlpqSxJqVM3jxncPUNztjsN7ZXkvAD5fPzo9xdBPLaFtQD+MM/+4B7rDWduJMml0P/FlkQhMRiY0Ur4c118+gKDc9tG/jzlp2HFCRmmga7Yq6HcA3AYwxWcaYHGvt9yMamYhIDCV7PU5ZpE1HOHHKGVC8aecJAgG4MlTDQCJp1Os8GWO+bow5DrQATcaYk8aYv4xcaCIisZXs9XDnigqm5J2p0/furhNsV7nPqBjtPKhHgO/jDDdfAVwP/CPwP4wx34xceCIiseVN8nDninIlqRgY7SCJbwB/aq39ddi+d40xR3Fq5/14zCMTEYkTQ0nqpY1VnDjl1AF4d9cJ3C6XBk5E0Gi7+PJxhniP9BHOpF0RkYR2tpbUxp21GoIeQaNNUHs4+/yiz+MUbxURSXhDSao4bHTfO9trNZk3QkbbxfcY8IoxZjmwJbhvOXArcG8kAhMRiUdOkqrgxY1HqGtyRvdt2F6D2+1ifkVujKNLLKNqQVlr1wKrgV6cWnz3AW3AUmvty5ELT0Qk/gwNQS8MWyp+/bYaKqtVYHYsjbYFhbX2HeCdCMYiIjJuDA1Bf/GdIzS0dAULzB7H7XZpqY4xcs4EFawM/iOc1lIv8ALwl2O1/pOIyHiXmpzEXSsq+N07hznV2k0gEODND47hdruYWZJz4ReQ8zpfF9/3gDuBf8BZUuMOnHp8IiISlBosMDu0VIc/EGDte0epPqnv8pfqfAnqPuBL1trvW2ufxBnFt8YY441OaCIi40NasMDs0KKH/kCAVzdXcby+PcaRjW/nS1AlDB9CvjV4fmFEIxIRGYd8qV7uXjmTrPRkAAb9AV55t4oTjVrg+2KdL0F5gMGhDWttAOdeVHKkgxIRGY8y0pwklZHmdDQNDPp5+d2q0LId8smMulisiIhcWFZ6MnevnIkv1UlSff2DvLjRGUQhn8yFhpl/2RgT3j5NAh4wxgybNm2t/fmYRyYiMk7lZKaw5voKXlh/mJ6+AXr7Bvn9O4e594aZTMpMjXV448b5EtQx4Osj9tUB/9eIfQFACUpEJExudlpoCHpf/yDdvQP8fsNhPnPjLDJ9ulMyGudMUNbasijGISKScAom+7jzugpefOcw/YN+Orr7+f0GpyU11AUo56Z7UCIiEVScl85t15ThcbsAaO3o5cWNR+jpG4hxZPFPCUpEJMKmFWVx81XTcbmcJHWqtZuXN1XRPzB4gSsnNiUoEZEomFGSw+olpaHtuqZOXt1czeCgP4ZRxTclKBGRKJlTNpkVV0wNbR+rb+eND47h9wdiGFX8UoISEYmiy2fls2x+UWj7UE0rG7bXEAgoSY2kBCUiEmVL5xZy+cz80PbeI028t6cuhhHFJyUoEZEoc7lcXHfFFOZMP7Nu1EeV9ew80BjDqOKPEpSISAy4XC5uWDKNsuKs0L6NO2upPKpVeYcoQYmIxIjH7eKWq8uYkpce2vfW1uMc1VpSgBKUiEhMeZPc3H5tObnZaUBwLakt1dQ1dcY2sDigBCUiEmNDS8cPrSU1MOjn5U1VNLf1xDiy2FKCEhGJA+lpXu5aMYO0FKdEak/fAC++c5iOrr4YRxY7SlAiInEiJzOFO6+rwJvkfDR3dPfz0gSu25cQCcoYc7kxZosxptMYs9sYszTWMYmIXIyCyT5uv6Ycd7C4bFNbD394t5qBCVgSadwnKGNMMvB74D+BHOBxYK0xJuu8F4qIxKnSwkw+tXRaaPvEqQ7eeP/ohCuJNO4TFLAK8Fprf2St7bfWPgfsBT4f27BERC7e7GmTuHbhlND24drTbNxRO6FKIiVCgpoH7B+xrxJYEINYRETGzJWmgCtmnymJtPvwKbbZhhhGFF2JkKAygK4R+7oAXwxiEREZU9cunMKs0jMlkbbsPjlhqk0kQoLqBNJG7PMBHTGIRURkTLlcLj61tJSSgozQvre2Hud4fXsMo4qOREhQ+wAzYt+c4H4RkXHP43Fz6/IycrNSgTPVJk61dsc2sAhLhAT1NuAyxnzLGOM1xnwBWAi8EOO4RETGTGpyEneuqCAjzQtAX/8gL286ktATecd9grLW9gG3AZ8BmoFvA3dba1W3XkQSSoYvmTtXVJDs9QCJP5E3KdYBjAVr7R7guljHISISabnZady2vIyXNh3B7w/Q1NbDa1uOcud15Xg8477NMUxi/TYiIhNAaWEmq5eUhrZrGtp5+6PEWzZeCUpEZBwy0ydz9WXFoe3Ko81s3V8fw4jGnhKUiMg4tXhOAfPKJ4e2P9hbR2V14syRUoISERmnXC4XKxeVUlqYGdr31kfHqWlIjDlSSlAiIuOYx+3ituVlZ1bk9TtzpFoSYLFDJSgRkXEu2evhzuvKSU915kj19g3y0qYjdPX0xziyS6MEJSKSADJ8ydxxbTne4FDzts4+/rB5fK8jpQQlIpIgCib7uPnq6bhczmKHdU2drNt6fNwOP1eCEhFJIOVTsrkubB2pg8db2LpvfA4/V4ISEUkwC2flcdmMvND2B/vqsONwiQ4lKBGRBONyubj+iqlMCx9+/uFxTpwaX6sQKUGJiCQgt9vFLcvLmBxcomPQH+DVzdWc7uiNcWSjpwQlIpKgUrwePn1dBWkpTl3w7t4BXnm3it7+wRhHNjpKUCIiCSwrPZnbrynH43ZG9jW39fD6e9X4/fE/sk8JSkQkwRXnpXNjWPXzY3XtvLvzRAwjGh0lKBGRCcBMn8ySuYWh7Z2HGtl7pCmGEV2YEpSIyARx1fwiZpTkhLY3bKuJ68KySlAiIhOEy+XiU0tLyc8JFpYNBHhty1Fa2+NzZJ8SlIjIBOJN8nDHteX4goVle/rid2SfEpSIyAST4Uvm9mvKQiP7Wtp7WPve0bgb2acEJSIyARXlDh/Zd7Sujff2nIxhRB+nBCUiMkGZ6ZNZZApC29tsQ1zV7FOCEhGZwK6+rJjy4qzQ9lsfHqe+uSuGEZ2hBCUiMoG53S5uumr6iJp9VXR2x341XiUoEZEJLtnrjOxLSfYA0NHdz6tbYr8arxKUiIiQnZHCrVeXDVuNd8O2mpiuxqsEJSIiAJQWZg5bjXd/dTO7Dp2KWTxKUCIiErJwVh5zyyaHtt/deSJm5ZCUoEREJMTlcrFyUQmFk33AmXJIbZ19UY9FCUpERIZJ8ri57Zrh5ZBe3VxF/0B0B00oQYmIyMdkpHm5bXkZ7mA5pMbWbt768HhUB00oQYmIyFkV56Wz8sqS0PbB4y3sONAYtfdPito7jZIx5lvASmvt3WH7pgG/BK4GGoD/21r7h+AxF/A/ga8BycC/AH9hrR2IduwiIolmfkUujS1d7Akubrh590nyctIoLcyM+HvHTQvKGJNhjHkS+OFZDj8H7AJygT8BnjPGVASPfQ24F1gEzAKWAn8V+YhFRCaGFVdMpSg3HYBAIMDr70Vn0ETcJCjgFaAc+EX4TmPMbGAJ8Ji1ts9a+xbwIvDV4CkPAT+y1tZYaxuBvwYejlrUIiIJzuNxc+vysqgPmohaF58xJhmYfJZDAWttPfBFa+0JY8xfA8Vhx+cBx6y1nWH7KoFlYcf3jTg2xRgz2VobP2V5RUTGsaFBEy9sOITfH6CxtZsN246zeum0UPWJsRbNFtQ1wMmzPGoBrLUnznFdBjCytG4X4DvH8aHnPkREZMwU56Wz4oqpoe3Koy3sOdwUsfeLWgvKWrseuJg02wmkjdjnAzrOcXwoMXUgIiJj6rKKXBqau9hf7XRQbdxRS15OGsV56WP+XvF0D+pc9gHTjDHhSWgOZ7r19gFmxLGT1trWKMUnIjJhDFWaKJh0ptLE+m01EXmvuE9Q1loL7AQeN8akGGNuANYAzwZP+RXwiDFmujEmD2eQxK9iEqyIyASQFBw0kZbidML19Q9GZAJv3M2DOofPAM/gzIE6BXzVWrsneOxpoBDYjNO99xvgsVgEKSIyUWSlJ/PZ1bM5cKyFsuKsiAyUcMVyrY94YIwpA6rWrVtHSUnJhU4XEZExVFNTw+rVqwHKrbXV4cfivotPREQmJiUoERGJS0pQIiISl5SgREQkLilBiYhIXFKCEhGRuDRe5kFFkgegrq4u1nGIiEw4YZ+9npHHlKCCldPvv//+WMchIjKRFQOHw3coQcFWYAVOZfXBGMciIjLReHCS09aRByZ8JQkREYlPGiQhIiJxSQlKRETikhKUiIjEJSUoERGJS0pQIiISl5SgREQkLilBiYhIXFKCEhGRuKRKEhfJGHM58DSwEDgCfMVa+7GZ0InKGHMT8H1gFtAAPGmt/UVso4o+Y0wOsAt4zFr7rzEOJ6qMMcXA/wJuAHqAZ6y134ltVNFjjLka+CfAAI3A9621/zu2UUWHMWYZ8LK1tiC4nQz8FLgPpyLPU9baJy71fdSCugjBf4zfA/8J5ACPA2uNMVkxDSxKjDGlwPPA3+L8/l8EnjDG3BLTwGLjaWBqrIOIkd/jlAgrBK4GHjLGfCm2IUWHMcaN8/v/k7U2G+dv4KfBL64JyxjjMsb8MbAWSA479D2cRD0DWIrzf+HBS30/JaiLswrwWmt/ZK3tt9Y+B+wFPh/bsKKmDHjWWvuCtdYfbDmuB66NaVRRZox5CMgCdsc6lmgzxlwFVAD/3VrbY62twvm7eDumgUXPJKAAcBljXEAAGAD6YhpV5H0P+DrOl9NwDwGPW2tbrLXVwA+Ahy/1zZSgLs48YP+IfZXAghjEEnXW2o3W2j8d2jbGTMYpuLs9dlFFlzGmHPgu8JVYxxIji3ES818bY2qNMYeBe6y1J2McV1RYa5twurT+P6Afp9DpX1lrR34uJJqnrbWLgQ+HdgS7uYuBfWHnjcnnoRLUxckAukbs6wJ8MYglpowx2cCLwPs4XR4JzxjjAf4deMRaO1EXEhv6UtKP05K6F3hkgnXx9QBfAtJwWo/fNcbcHMu4Is1ae+IsuzOCP8M/E8fk81CDJC5OJ85/ynA+oCMGscSMMWY2TlLaB9xvrfXHOKRo+Q5grbW/jXUgMdQLtFlr/zq4vdMY879xEtWzMYsqeu4FrrXW/kVwe4Mx5pc43VprYxdWTHQGf4Z/Jo7J56FaUBdnH84NwXBzGN7ETWjGmOtxWk2/A+6z1vbEOKRo+gJwnzGm1RjTitOV8XNjzM9jHFc0VQK+4IChIRPpC28pkDJi3wBOi3JCsda2AHUM/0wck8/DifQfaiy9jXNz9Fs4/dCfwRlu/kJMo4oSY8wM4GXg29ban8Q6nmiz1s4J3zbG7AB+NMGGmb+BM7T6h8aYP8f5cPoqzg30iWAtzsjVrwH/DCwC/gT445hGFTu/wuni3IXT5fcI8ONLfVG1oC6CtbYPuA0nMTUD3wbuttY2xjSw6PkGkInzB9oR9vj7WAcm0RFsMa/Euf90EngN+Adr7fMxDSxKrLV7cbr5HgZacbo1/9JaOyHuw57FY8AenNHMW3GmoTx9qS+qFXVFRCQuqQUlIiJxSQlKRETikhKUiIjEJSUoERGJS0pQIiISl5SgREQkLmmirkiEGGP+FafK87l8D6cK/NtAprU2KqWygrUE3wUetNYeOM95buA94I+stTYasYmEUwtKJHK+iVPluRinmCjAsrB9PwA2B593nuX6SPnvwM7zJSeAYG3Fv2EMJlyKXAxN1BWJAmPMZTjLU5QH18uJVRypwDHgRmvtnlFecxj4qrV2fSRjExlJXXwiMWSMWUVYF58xJoCzOuujOPXtPgQeAP4C+COgDXjUWvur4PWZwA9xltoOAG8B3zzHsgjgFLptDU9OxpjvAF8D8nHWOfsra+2rYde8gNMaXD8Gv7LIqKmLTyT+fB/4M5xl1KcB23AS01Lgt8AvjDFDa/A8g5PIbsGpjRcAXjfGnOvL5x04dfMAMMbcE3yvB3AqUL8C/MYYkxV2zWvAp87zmiIRoQQlEn9+Zq1921q7A6dqfAdOq8YCT+Gsu1NujKnAaRF9yVq7Ndgq+iOgDLj1HK+9BKeg55AynLWdjga7Hv8Gpwhq+LIR+3AqVA+r4i4SafpGJBJ/DoU97wKqrbVDN4uH1t1KAaYHn1tjhi1P5sNpVb18ltcuBE6Fbf87zkjDI8aYj3BWR/4Xa2132DlNwZ8Fn/D3ELkkakGJxJ+Ri96da6XipOC5VwJXhD1mA/9yjmv8gGtoI7hEzGKcFtdm4MvAruCgjiFDnxODo/4NRMaAEpTI+LUf8ALp1tpD1tpDOGszPYmTpM6mDmcwBADGmHuBh621a62138RpebUDt4ddkx92rUjUqItPZJyy1lpjzIvAvxljvoGzwu3jOIMrKs9x2UfA5WHbHuBJY0w9zojBq4Gi4PMhlwMtDO96FIk4taBExreHcJLJ73BWMs0GbrLWtp7j/FdwRvsBYK39DfBdnFbXAeBvgf9mrX0r7Jrrgdesterik6jSRF2RCcQY4wOqgVuttdtGcb4bOIozUnBjhMMTGUYtKJEJxFrbhdNa+sYoL1kDHFFyklhQghKZeP4RWGhGjE0fKdh6+jbwp1GJSmQEdfGJiEhcUgtKRETikhKUiIjEJSUoERGJS0pQIiISl5SgREQkLv3/oFvnnUA0aHgAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "def plot_position(results):\n",
    "    plot(results.y, label='y')\n",
    "    decorate(xlabel='Time (s)',\n",
    "             ylabel='Position (m)')\n",
    "\n",
    "plot_position(results)\n",
    "savefig('figs/chap09-fig01.pdf')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Onto the sidewalk\n",
    "\n",
    "To figure out when the penny hit the sidewalk, we can use `crossings`, which finds the times where a `Series` passes through a given value."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "metadata": {},
   "outputs": [],
   "source": [
    "def crossings(series, value):\n",
    "    \"\"\"Find the labels where the series passes through value.\n",
    "\n",
    "    The labels in series must be increasing numerical values.\n",
    "\n",
    "    series: Series\n",
    "    value: number\n",
    "\n",
    "    returns: sequence of labels\n",
    "    \"\"\"\n",
    "    units = get_units(series.values[0])\n",
    "    values = magnitudes(series - value)\n",
    "    interp = InterpolatedUnivariateSpline(series.index, values)\n",
    "    return interp.roots()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([8.81788535])"
      ]
     },
     "execution_count": 35,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "t_crossings = crossings(results.y, 0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For this example there should be just one crossing, the time when the penny hits the sidewalk."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 36,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "8.81788534972054 second"
      ],
      "text/latex": [
       "$8.81788534972054\\ \\mathrm{second}$"
      ],
      "text/plain": [
       "8.81788534972054 <Unit('second')>"
      ]
     },
     "execution_count": 36,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "t_sidewalk = t_crossings[0] * s"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can compare that to the exact result.  Without air resistance, we have\n",
    "\n",
    "$v = -g t$\n",
    "\n",
    "and\n",
    "\n",
    "$y = 381 - g t^2 / 2$\n",
    "\n",
    "Setting $y=0$ and solving for $t$ yields\n",
    "\n",
    "$t = \\sqrt{\\frac{2 y_{init}}{g}}$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "8.817885349720552 second"
      ],
      "text/latex": [
       "$8.817885349720552\\ \\mathrm{second}$"
      ],
      "text/plain": [
       "8.817885349720552 <Unit('second')>"
      ]
     },
     "execution_count": 37,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "sqrt(2 * init.y / g)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The estimate is accurate to about 10 decimal places."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Events\n",
    "\n",
    "Instead of running the simulation until the penny goes through the sidewalk, it would be better to detect the point where the penny hits the sidewalk and stop.  `run_ralston` provides exactly the tool we need, **event functions**.\n",
    "\n",
    "Here's an event function that returns the height of the penny above the sidewalk:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 38,
   "metadata": {},
   "outputs": [],
   "source": [
    "def event_func(state, t, system):\n",
    "    \"\"\"Return the height of the penny above the sidewalk.\n",
    "    \"\"\"\n",
    "    y, v = state\n",
    "    return y"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And here's how we pass it to `run_ralston`.  The solver should run until the event function returns 0, and then terminate."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 39,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>values</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>message</th>\n",
       "      <td>Success</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "message    Success\n",
       "dtype: object"
      ]
     },
     "execution_count": 39,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "results, details = run_ralston(system, slope_func, events=event_func)\n",
    "details"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The message from the solver indicates the solver stopped because the event we wanted to detect happened.\n",
    "\n",
    "Here are the results:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 40,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>y</th>\n",
       "      <th>v</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0.000000</th>\n",
       "      <td>381 meter</td>\n",
       "      <td>0.0 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>0.100000</th>\n",
       "      <td>380.951 meter</td>\n",
       "      <td>-0.9800000000000001 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>0.200000</th>\n",
       "      <td>380.80400000000003 meter</td>\n",
       "      <td>-1.9600000000000002 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>0.300000</th>\n",
       "      <td>380.559 meter</td>\n",
       "      <td>-2.9400000000000004 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>0.400000</th>\n",
       "      <td>380.216 meter</td>\n",
       "      <td>-3.9200000000000004 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>0.500000</th>\n",
       "      <td>379.77500000000003 meter</td>\n",
       "      <td>-4.9 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>0.600000</th>\n",
       "      <td>379.23600000000005 meter</td>\n",
       "      <td>-5.880000000000001 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>0.700000</th>\n",
       "      <td>378.59900000000005 meter</td>\n",
       "      <td>-6.860000000000001 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>0.800000</th>\n",
       "      <td>377.86400000000003 meter</td>\n",
       "      <td>-7.840000000000002 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>0.900000</th>\n",
       "      <td>377.031 meter</td>\n",
       "      <td>-8.820000000000002 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1.000000</th>\n",
       "      <td>376.1 meter</td>\n",
       "      <td>-9.800000000000002 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1.100000</th>\n",
       "      <td>375.071 meter</td>\n",
       "      <td>-10.780000000000003 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1.200000</th>\n",
       "      <td>373.944 meter</td>\n",
       "      <td>-11.760000000000003 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1.300000</th>\n",
       "      <td>372.719 meter</td>\n",
       "      <td>-12.740000000000004 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1.400000</th>\n",
       "      <td>371.396 meter</td>\n",
       "      <td>-13.720000000000004 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1.500000</th>\n",
       "      <td>369.975 meter</td>\n",
       "      <td>-14.700000000000005 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1.600000</th>\n",
       "      <td>368.456 meter</td>\n",
       "      <td>-15.680000000000005 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1.700000</th>\n",
       "      <td>366.839 meter</td>\n",
       "      <td>-16.660000000000004 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1.800000</th>\n",
       "      <td>365.124 meter</td>\n",
       "      <td>-17.640000000000004 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1.900000</th>\n",
       "      <td>363.31100000000004 meter</td>\n",
       "      <td>-18.620000000000005 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2.000000</th>\n",
       "      <td>361.40000000000003 meter</td>\n",
       "      <td>-19.600000000000005 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2.100000</th>\n",
       "      <td>359.391 meter</td>\n",
       "      <td>-20.580000000000005 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2.200000</th>\n",
       "      <td>357.284 meter</td>\n",
       "      <td>-21.560000000000006 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2.300000</th>\n",
       "      <td>355.079 meter</td>\n",
       "      <td>-22.540000000000006 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2.400000</th>\n",
       "      <td>352.776 meter</td>\n",
       "      <td>-23.520000000000007 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2.500000</th>\n",
       "      <td>350.375 meter</td>\n",
       "      <td>-24.500000000000007 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2.600000</th>\n",
       "      <td>347.876 meter</td>\n",
       "      <td>-25.480000000000008 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2.700000</th>\n",
       "      <td>345.279 meter</td>\n",
       "      <td>-26.460000000000008 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2.800000</th>\n",
       "      <td>342.584 meter</td>\n",
       "      <td>-27.44000000000001 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2.900000</th>\n",
       "      <td>339.791 meter</td>\n",
       "      <td>-28.42000000000001 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>...</th>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6.000000</th>\n",
       "      <td>204.60000000000005 meter</td>\n",
       "      <td>-58.799999999999926 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6.100000</th>\n",
       "      <td>198.67100000000005 meter</td>\n",
       "      <td>-59.77999999999992 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6.200000</th>\n",
       "      <td>192.64400000000006 meter</td>\n",
       "      <td>-60.75999999999992 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6.300000</th>\n",
       "      <td>186.51900000000006 meter</td>\n",
       "      <td>-61.73999999999992 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6.400000</th>\n",
       "      <td>180.29600000000008 meter</td>\n",
       "      <td>-62.719999999999914 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6.500000</th>\n",
       "      <td>173.97500000000008 meter</td>\n",
       "      <td>-63.69999999999991 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6.600000</th>\n",
       "      <td>167.5560000000001 meter</td>\n",
       "      <td>-64.67999999999991 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6.700000</th>\n",
       "      <td>161.0390000000001 meter</td>\n",
       "      <td>-65.65999999999991 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6.800000</th>\n",
       "      <td>154.42400000000012 meter</td>\n",
       "      <td>-66.63999999999992 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6.900000</th>\n",
       "      <td>147.71100000000013 meter</td>\n",
       "      <td>-67.61999999999992 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7.000000</th>\n",
       "      <td>140.90000000000015 meter</td>\n",
       "      <td>-68.59999999999992 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7.100000</th>\n",
       "      <td>133.99100000000016 meter</td>\n",
       "      <td>-69.57999999999993 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7.200000</th>\n",
       "      <td>126.98400000000017 meter</td>\n",
       "      <td>-70.55999999999993 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7.300000</th>\n",
       "      <td>119.87900000000018 meter</td>\n",
       "      <td>-71.53999999999994 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7.400000</th>\n",
       "      <td>112.67600000000019 meter</td>\n",
       "      <td>-72.51999999999994 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7.500000</th>\n",
       "      <td>105.3750000000002 meter</td>\n",
       "      <td>-73.49999999999994 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7.600000</th>\n",
       "      <td>97.9760000000002 meter</td>\n",
       "      <td>-74.47999999999995 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7.700000</th>\n",
       "      <td>90.4790000000002 meter</td>\n",
       "      <td>-75.45999999999995 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7.800000</th>\n",
       "      <td>82.8840000000002 meter</td>\n",
       "      <td>-76.43999999999996 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7.900000</th>\n",
       "      <td>75.1910000000002 meter</td>\n",
       "      <td>-77.41999999999996 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8.000000</th>\n",
       "      <td>67.4000000000002 meter</td>\n",
       "      <td>-78.39999999999996 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8.100000</th>\n",
       "      <td>59.51100000000021 meter</td>\n",
       "      <td>-79.37999999999997 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8.200000</th>\n",
       "      <td>51.524000000000214 meter</td>\n",
       "      <td>-80.35999999999997 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8.300000</th>\n",
       "      <td>43.43900000000022 meter</td>\n",
       "      <td>-81.33999999999997 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8.400000</th>\n",
       "      <td>35.25600000000022 meter</td>\n",
       "      <td>-82.31999999999998 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8.500000</th>\n",
       "      <td>26.97500000000022 meter</td>\n",
       "      <td>-83.29999999999998 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8.600000</th>\n",
       "      <td>18.596000000000224 meter</td>\n",
       "      <td>-84.27999999999999 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8.700000</th>\n",
       "      <td>10.119000000000225 meter</td>\n",
       "      <td>-85.25999999999999 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8.800000</th>\n",
       "      <td>1.5440000000002243 meter</td>\n",
       "      <td>-86.24 meter / second</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8.817802</th>\n",
       "      <td>-4.440892098500626e-16 meter</td>\n",
       "      <td>-86.41446327683617 meter / second</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>90 rows × 2 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "                                     y                                   v\n",
       "0.000000                     381 meter                  0.0 meter / second\n",
       "0.100000                 380.951 meter  -0.9800000000000001 meter / second\n",
       "0.200000      380.80400000000003 meter  -1.9600000000000002 meter / second\n",
       "0.300000                 380.559 meter  -2.9400000000000004 meter / second\n",
       "0.400000                 380.216 meter  -3.9200000000000004 meter / second\n",
       "0.500000      379.77500000000003 meter                 -4.9 meter / second\n",
       "0.600000      379.23600000000005 meter   -5.880000000000001 meter / second\n",
       "0.700000      378.59900000000005 meter   -6.860000000000001 meter / second\n",
       "0.800000      377.86400000000003 meter   -7.840000000000002 meter / second\n",
       "0.900000                 377.031 meter   -8.820000000000002 meter / second\n",
       "1.000000                   376.1 meter   -9.800000000000002 meter / second\n",
       "1.100000                 375.071 meter  -10.780000000000003 meter / second\n",
       "1.200000                 373.944 meter  -11.760000000000003 meter / second\n",
       "1.300000                 372.719 meter  -12.740000000000004 meter / second\n",
       "1.400000                 371.396 meter  -13.720000000000004 meter / second\n",
       "1.500000                 369.975 meter  -14.700000000000005 meter / second\n",
       "1.600000                 368.456 meter  -15.680000000000005 meter / second\n",
       "1.700000                 366.839 meter  -16.660000000000004 meter / second\n",
       "1.800000                 365.124 meter  -17.640000000000004 meter / second\n",
       "1.900000      363.31100000000004 meter  -18.620000000000005 meter / second\n",
       "2.000000      361.40000000000003 meter  -19.600000000000005 meter / second\n",
       "2.100000                 359.391 meter  -20.580000000000005 meter / second\n",
       "2.200000                 357.284 meter  -21.560000000000006 meter / second\n",
       "2.300000                 355.079 meter  -22.540000000000006 meter / second\n",
       "2.400000                 352.776 meter  -23.520000000000007 meter / second\n",
       "2.500000                 350.375 meter  -24.500000000000007 meter / second\n",
       "2.600000                 347.876 meter  -25.480000000000008 meter / second\n",
       "2.700000                 345.279 meter  -26.460000000000008 meter / second\n",
       "2.800000                 342.584 meter   -27.44000000000001 meter / second\n",
       "2.900000                 339.791 meter   -28.42000000000001 meter / second\n",
       "...                                ...                                 ...\n",
       "6.000000      204.60000000000005 meter  -58.799999999999926 meter / second\n",
       "6.100000      198.67100000000005 meter   -59.77999999999992 meter / second\n",
       "6.200000      192.64400000000006 meter   -60.75999999999992 meter / second\n",
       "6.300000      186.51900000000006 meter   -61.73999999999992 meter / second\n",
       "6.400000      180.29600000000008 meter  -62.719999999999914 meter / second\n",
       "6.500000      173.97500000000008 meter   -63.69999999999991 meter / second\n",
       "6.600000       167.5560000000001 meter   -64.67999999999991 meter / second\n",
       "6.700000       161.0390000000001 meter   -65.65999999999991 meter / second\n",
       "6.800000      154.42400000000012 meter   -66.63999999999992 meter / second\n",
       "6.900000      147.71100000000013 meter   -67.61999999999992 meter / second\n",
       "7.000000      140.90000000000015 meter   -68.59999999999992 meter / second\n",
       "7.100000      133.99100000000016 meter   -69.57999999999993 meter / second\n",
       "7.200000      126.98400000000017 meter   -70.55999999999993 meter / second\n",
       "7.300000      119.87900000000018 meter   -71.53999999999994 meter / second\n",
       "7.400000      112.67600000000019 meter   -72.51999999999994 meter / second\n",
       "7.500000       105.3750000000002 meter   -73.49999999999994 meter / second\n",
       "7.600000        97.9760000000002 meter   -74.47999999999995 meter / second\n",
       "7.700000        90.4790000000002 meter   -75.45999999999995 meter / second\n",
       "7.800000        82.8840000000002 meter   -76.43999999999996 meter / second\n",
       "7.900000        75.1910000000002 meter   -77.41999999999996 meter / second\n",
       "8.000000        67.4000000000002 meter   -78.39999999999996 meter / second\n",
       "8.100000       59.51100000000021 meter   -79.37999999999997 meter / second\n",
       "8.200000      51.524000000000214 meter   -80.35999999999997 meter / second\n",
       "8.300000       43.43900000000022 meter   -81.33999999999997 meter / second\n",
       "8.400000       35.25600000000022 meter   -82.31999999999998 meter / second\n",
       "8.500000       26.97500000000022 meter   -83.29999999999998 meter / second\n",
       "8.600000      18.596000000000224 meter   -84.27999999999999 meter / second\n",
       "8.700000      10.119000000000225 meter   -85.25999999999999 meter / second\n",
       "8.800000      1.5440000000002243 meter               -86.24 meter / second\n",
       "8.817802  -4.440892098500626e-16 meter   -86.41446327683617 meter / second\n",
       "\n",
       "[90 rows x 2 columns]"
      ]
     },
     "execution_count": 40,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "results"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "With the `events` option, the solver returns the actual time steps it computed, which are not necessarily equally spaced. \n",
    "\n",
    "The last time step is when the event occurred:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 41,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "8.81780237518735 second"
      ],
      "text/latex": [
       "$8.81780237518735\\ \\mathrm{second}$"
      ],
      "text/plain": [
       "8.81780237518735 <Unit('second')>"
      ]
     },
     "execution_count": 41,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "t_sidewalk = get_last_label(results) * s"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The result is accurate to about 15 decimal places.\n",
    "\n",
    "We can also check the velocity of the penny when it hits the sidewalk:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 42,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "-86.41446327683617 meter/second"
      ],
      "text/latex": [
       "$-86.41446327683617\\ \\frac{\\mathrm{meter}}{\\mathrm{second}}$"
      ],
      "text/plain": [
       "-86.41446327683617 <Unit('meter / second')>"
      ]
     },
     "execution_count": 42,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "v_sidewalk = get_last_value(results.v)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And convert to kilometers per hour."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 43,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "-311.09206779661025 kilometer/hour"
      ],
      "text/latex": [
       "$-311.09206779661025\\ \\frac{\\mathrm{kilometer}}{\\mathrm{hour}}$"
      ],
      "text/plain": [
       "-311.09206779661025 <Unit('kilometer / hour')>"
      ]
     },
     "execution_count": 43,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "km = UNITS.kilometer\n",
    "h = UNITS.hour\n",
    "v_sidewalk.to(km / h)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "If there were no air resistance, the penny would hit the sidewalk (or someone's head) at more than 300 km/h.\n",
    "\n",
    "So it's a good thing there is air resistance."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
